chi2_shift¶
- image_registration.chi2_shifts.chi2_shift(im1, im2, err=None, upsample_factor='auto', boundary='wrap', nthreads=1, use_numpy_fft=False, zeromean=False, nfitted=2, verbose=False, return_error=True, return_chi2array=False, max_auto_size=512, max_nsig=1.1)[source]¶
Find the offsets between image 1 and image 2 using the DFT upsampling method (http://www.mathworks.com/matlabcentral/fileexchange/18401-efficient-subpixel-image-registration-by-cross-correlation/content/html/efficient_subpixel_registration.html) combined with \(\chi^2\) to measure the errors on the fit
Equation 1 gives the \(\chi^2\) value as a function of shift, where Y is the model as a function of shift:
\[\chi^2(dx,dy) = \Sigma_{ij} \frac{(X_{ij}-Y_{ij}(dx,dy))^2}{\sigma_{ij}^2}\]Equation 2-4:
\begin{align} \mathrm{Term~1:} & f(dx,dy) & = & \Sigma_{ij} \frac{X_{ij}^2}{\sigma_{ij}^2} \\ & f(dx,dy) & = & f(0,0) , \forall dx,dy \\ \mathrm{Term~2:} & g(dx,dy) & = & -2 \Sigma_{ij} \frac{X_{ij}Y_{ij}(dx,dy)}{\sigma_{ij}^2} = -2 \Sigma_{ij} \left(\frac{X_{ij}}{\sigma_{ij}^2}\right) Y_{ij}(dx,dy) \\ \mathrm{Term~3:} & h(dx,dy) & = & \Sigma_{ij} \frac{Y_{ij}(dx,dy)^2}{\sigma_{ij}^2} = \Sigma_{ij} \left(\frac{1}{\sigma_{ij}^2}\right) Y^2_{ij}(dx,dy) \end{align}The cross-correlation can be computed with fourier transforms, and is defined
\[CC_{m,n}(x,y) = \Sigma_{ij} x^*_{ij} y_{(n+i)(m+j)}\]which can then be applied to our problem, noting that the cross-correlation has the same form as term 2 and 3 in \(\chi^2\) (term 1 is a constant, with no dependence on the shift)
\begin{align} \mathrm{Term~2:} & CC(X/\sigma^2,Y)[dx,dy] & = & \Sigma_{ij} \left(\frac{X_{ij}}{\sigma_{ij}^2}\right)^* Y_{ij}(dx,dy) \\ \mathrm{Term~3:} & CC(\sigma^{-2},Y^2)[dx,dy] & = & \Sigma_{ij} \left(\frac{1}{\sigma_{ij}^2}\right)^* Y^2_{ij}(dx,dy) \end{align}Technically, only terms 2 and 3 has any effect on the resulting image, since term 1 is the same for all shifts, and the quantity of interest is \(\Delta \chi^2\) when determining the best-fit shift and error.
The resulting shifts are limited to an accuracy of +/-0.5 pixels in the upsampled image frame. That is not a Gaussian uncertainty but a quantized limit: the best solution will always be +/-0.5 pixels offset because we’re zooming in on an even grid, so the “best” position is required to be a discrete pixel center. If you’re looking at an image with exactly zero shift, it will have exactly +/- 1/usfac/2 systematic error in the resulting solution.
- Parameters:
- im1np.ndarray
- im2np.ndarray
The images to register.
- errnp.ndarray
Per-pixel error in image 2
- boundary‘wrap’,’constant’,’reflect’,’nearest’
Option to pass to map_coordinates for determining what to do with shifts outside of the boundaries.
- upsample_factorint or ‘auto’
upsampling factor; governs accuracy of fit (1/usfac is best accuracy) (can be “automatically” determined based on chi^2 error)
- return_errorbool
Returns the “fit error” (1-sigma in x and y) based on the delta-chi2 values
- return_chi2_arraybool
Returns the x and y shifts and the chi2 as a function of those shifts in addition to other returned parameters. i.e., the last return from this function will be a tuple (x, y, chi2)
- zeromeanbool
Subtract the mean from the images before cross-correlating? If no, you may get a 0,0 offset because the DC levels are strongly correlated.
- verbosebool
Print error message if upsampling factor is inadequate to measure errors
- use_numpy_fftbool
Force use numpy’s fft over fftw? (only matters if you have fftw installed)
- nthreadsbool
Number of threads to use for fft (only matters if you have fftw installed)
- nfittedint
number of degrees of freedom in the fit (used for chi^2 computations). Should probably always be 2.
- max_auto_sizeint
Maximum zoom image size to create when using auto-upsampling
- Returns:
- dx,dyfloat,float
Measures the amount im2 is offset from im1 (i.e., shift im2 by -1 * these #’s to match im1)
- errx,erryfloat,float
optional, error in x and y directions
- xvals,yvals,chi2n_upsampledndarray,ndarray,ndarray,
x,y positions (in original chi^2 coordinates) of the chi^2 values and their corresponding chi^2 value
Examples
Create a 2d array, shift it in both directions, then use chi2_shift to determine the shift
>>> import image_registration >>> rr = ((np.indices([100,100]) - np.array([50.,50.])[:,None,None])**2).sum(axis=0)**0.5 >>> image = np.exp(-rr**2/(3.**2*2.)) * 20 >>> shifted = np.roll(np.roll(image,12,0),5,1) + np.random.randn(100,100) >>> dx,dy,edx,edy = chi2_shift(image, shifted, upsample_factor='auto') >>> shifted2 = image_registration.fft_tools.shift2d(image,3.665,-4.25) + np.random.randn(100,100) >>> dx2,dy2,edx2,edy2 = chi2_shift(image, shifted2, upsample_factor='auto')