In [1]:

    

import numpy as np
import image_registration
%load_ext autoreload
%autoreload 2
%aimport image_registration


    

    




WARNING: Module py was already imported from /Users/adam/virtual-python/lib/python2.7/site-packages/astropy/extern/pytest.pyc/py, but /Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages is being added to sys.path [pkg_resources]
WARNING: Module _pytest was already imported from /Users/adam/virtual-python/lib/python2.7/site-packages/astropy/extern/pytest.pyc/_pytest, but /Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages is being added to sys.path [pkg_resources]
WARNING: Module pytest was already imported from /Users/adam/virtual-python/lib/python2.7/site-packages/astropy/extern/pytest.pyc/pytest, but /Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages is being added to sys.path [pkg_resources]

In [2]:

    

im = image_registration.tests.make_extended(100)


    

    




WARNING: RuntimeWarning: divide by zero encountered in power [image_registration.tests.registration_testing]
WARNING: RuntimeWarning: divide by zero encountered in power [image_registration.tests.registration_testing]
WARNING: RuntimeWarning: invalid value encountered in multiply [image_registration.tests.registration_testing]

In [3]:

    

yy,xx = np.indices(im.shape)


    

In [4]:

    

figure()
imshow(xx)
figure()
imshow(yy)
figure()
imshow(im)


    

    
Out[4]:





<matplotlib.image.AxesImage at 0x107f7a910>






    













    













    

In [5]:

    

im2 = image_registration.tests.make_offset_extended(im, -11.7, 10.9)


    

In [6]:

    

subplot(121);imshow(im)
subplot(122);imshow(im2)


    

    
Out[6]:





<matplotlib.image.AxesImage at 0x107f80ad0>






    

In [7]:

    

xc = image_registration.fft_tools.correlate2d(im,im2)


    

In [8]:

    

imshow(xc)


    

    
Out[8]:





<matplotlib.image.AxesImage at 0x108217d10>






    

In [8]:

    




    

In [8]:

    




    

In [9]:

    

i1,i2,t = image_registration.tests.make_offset_images(5,-7,100)


    

In [10]:

    

subplot(121);imshow(i1)
subplot(122);imshow(i2)


    

    
Out[10]:





<matplotlib.image.AxesImage at 0x1083fae50>






    

In [11]:

    

import astropy.nddata


    

In [12]:

    

astropy.nddata.convolve_fft


    

    
Out[12]:





<function astropy.nddata.convolution.convolve.convolve_fft>

In [14]:

    

#image_registration.tests.registration_testing.test_extended_shifts(1,1,100,False)


    

In [15]:

    

image_registration.tests.registration_testing.fit_extended_shifts(1,1,100,True)


    

    




Cannot estimate errors: need higher upsample factor.






    
Out[15]:





(-46.0, 37.00390625, 0.00390625, 0.00390625, 4)

In [22]:

    

reload(image_registration.tests.registration_testing)
reload(image_registration.tests)
reload(image_registration)
reload(image_registration.chi2_shifts)


    

    
Out[22]:





<module 'image_registration.chi2_shifts' from '/Users/adam/virtual-python/lib/python2.7/site-packages/image_registration/chi2_shifts.pyc'>

In [23]:

    

noise_taper=True


    

In [24]:

    

imsize=99
xsh=1
ysh=1
image = image_registration.tests.make_extended(imsize)
offset_image = image_registration.tests.make_offset_extended(image, xsh, ysh, noise_taper=noise_taper)
if noise_taper:
    noise = 1.0/image_registration.tests.edge_weight(imsize)
else:
    noise = 1.0
xoff,yoff,exoff,eyoff,(x,y,c2a) = image_registration.chi2_shift(image,offset_image,1.,return_error=True,verbose=True,upsample_factor='auto',return_chi2array=True)
print "SCALAR error: ",xoff,yoff,exoff,eyoff
print
xoff,yoff,exoff,eyoff,(x,y,c2) = image_registration.chi2_shift(image,offset_image,noise,return_error=True,verbose=True,upsample_factor='auto',return_chi2array=True)
c2map,term1,term2,term3 = image_registration.chi2n_map(image,offset_image,noise,return_all=True)

c2mapA,term1A,term2A,term3A = image_registration.chi2n_map(image,offset_image,1.,return_all=True)
print xoff,yoff,exoff,eyoff


    

    




Coarse xmax/ymax = 50,50, for offset 1.000000,1.000000
Minimum chi2: 18664.3   Max delta-chi2 (lowres): 45699.3  Min delta-chi2: 363.014
Selected upsample factor 256.0 for image size 512 and zoom factor 2.0 (max-sigma range was 1 for area 1)
Minimum chi2_ups: -2.13701e+09   Max delta-chi2: 3.93845e+07  Min delta-chi2: 23673.9
Cannot estimate errors: need higher upsample factor.
SCALAR error:  -0 -0 0.00390625 0.00390625

Coarse xmax/ymax = 64,19, for offset 15.000000,-30.000000
Minimum chi2: 5.19103e+09   Max delta-chi2 (lowres): 3.93778e+09  Min delta-chi2: 2.40917e+06
Selected upsample factor 256.0 for image size 512 and zoom factor 2.0 (max-sigma range was 1 for area 1)
Minimum chi2_ups: -1.28196e+09   Max delta-chi2: 8.30733e+06  Min delta-chi2: 7283.39
Cannot estimate errors: need higher upsample factor.
-15.99609375 31.0 0.00390625 0.00390625

In [19]:

    

subplot(131); imshow(image)
subplot(132); imshow(offset_image)
subplot(133); imshow(noise)
figure()
imshow(c2.real); colorbar()
figure()
imshow(term2); colorbar()
figure()
imshow(term3); colorbar()
figure()
imshow(c2map); colorbar()
figure()
subplot(131); imshow(image/noise**2); colorbar()
subplot(132); imshow(offset_image**2); colorbar()
subplot(133); imshow(noise**2); colorbar()


    

    
Out[19]:





<matplotlib.colorbar.Colorbar instance at 0x10891b7e8>






    













    













    













    













    













    

In [20]:

    

figure()
imshow(c2a.real); colorbar()
figure()
imshow(term2A); colorbar()
figure()
print term1A,term3A
imshow(c2mapA); colorbar()


    

    




121507.594882 96720.0720589






    
Out[20]:





<matplotlib.colorbar.Colorbar instance at 0x10900b878>






    













    













    

In [20]:

    




    

In [20]:

    




    

In [25]:

    

x = linspace(-5,5)


    

In [28]:

    

g = exp(-x**2)


    

In [29]:

    

plot(x,g)


    

    
Out[29]:





[<matplotlib.lines.Line2D at 0x1086d3b90>]






    

In [30]:

    

d = np.ones([3,3])


    

In [31]:

    

fft2(d)


    

    
Out[31]:





array([[ 9.+0.j,  0.+0.j,  0.+0.j],
       [ 0.+0.j,  0.+0.j,  0.+0.j],
       [ 0.+0.j,  0.+0.j,  0.+0.j]])

In [33]:

    

fft.fft(d,axis=0)


    

    
Out[33]:





array([[ 3.+0.j,  3.+0.j,  3.+0.j],
       [ 0.+0.j,  0.+0.j,  0.+0.j],
       [ 0.+0.j,  0.+0.j,  0.+0.j]])

In [34]:

    

fft.fft(fft.fft(d,axis=0),axis=1)


    

    
Out[34]:





array([[ 9.+0.j,  0.+0.j,  0.+0.j],
       [ 0.+0.j,  0.+0.j,  0.+0.j],
       [ 0.+0.j,  0.+0.j,  0.+0.j]])

In [1]:

    

x = linspace(-5,5,10)
g = exp(-x**2)


    

In [2]:

    

plot(g)


    

    
Out[2]:





[<matplotlib.lines.Line2D at 0x105892fd0>]






    

In [4]:

    

from numpy.fft import *


    

In [5]:

    

gft = fft(g)


    

In [10]:

    

guft = np.zeros(g.size*10,dtype=complex)
guft[:5] = gft[:5]
guft[-5:] = gft[-5:]


    

In [11]:

    

plot(guft.real)
plot(guft.imag)


    

    
Out[11]:





[<matplotlib.lines.Line2D at 0x105cfa050>]






    

In [25]:

    

gus=real(ifft(guft))*g.size
plot(linspace(x.min(),x.max(),guft.size),gus)
plot(x,g)


    

    
Out[25]:





[<matplotlib.lines.Line2D at 0x106fa68d0>]






    

In [18]:

    

gft,gft[-5:],gft[:5]


    

    
Out[18]:





(array([  1.59413217e+00 +0.00000000e+00j,
        -1.39813625e+00 -4.54282005e-01j,
         9.29589826e-01 +6.75386542e-01j,
        -4.37972931e-01 -6.02818024e-01j,
         1.09453270e-01 +3.36862526e-01j,
        -3.60822483e-16 -9.43689571e-16j,
         1.09453270e-01 -3.36862526e-01j,
        -4.37972931e-01 +6.02818024e-01j,
         9.29589826e-01 -6.75386542e-01j,  -1.39813625e+00 +4.54282005e-01j]),
 array([ -3.60822483e-16 -9.43689571e-16j,
         1.09453270e-01 -3.36862526e-01j,
        -4.37972931e-01 +6.02818024e-01j,
         9.29589826e-01 -6.75386542e-01j,  -1.39813625e+00 +4.54282005e-01j]),
 array([ 1.59413217+0.j        , -1.39813625-0.45428201j,
        0.92958983+0.67538654j, -0.43797293-0.60281802j,
        0.10945327+0.33686253j]))

In [21]:

    

x2 = linspace(-5,5,11)
g2 = exp(-x2**2)
plot(g2)


    

    
Out[21]:





[<matplotlib.lines.Line2D at 0x106cdcc50>]






    

In [22]:

    

gft2 = fft(g2)
guft2 = np.zeros(g2.size*10,dtype=complex)
guft2[:5] = gft2[:5]
guft2[-5:] = gft2[-5:]


    

In [23]:

    

plot(guft2.real)
plot(guft2.imag)


    

    
Out[23]:





[<matplotlib.lines.Line2D at 0x10677db90>]






    

In [26]:

    

gus2 = real(ifft(guft2))*g2.size
plot(linspace(x2.min(),x2.max(),guft2.size),gus2)
plot(x2,g2)


    

    
Out[26]:





[<matplotlib.lines.Line2D at 0x106feff90>]






    

In [28]:

    

plot(g2)
plot(gus2[::10])


    

    
Out[28]:





[<matplotlib.lines.Line2D at 0x106ff2290>]






    

In [29]:

    

plot(g)
plot(gus[::10])


    

    
Out[29]:





[<matplotlib.lines.Line2D at 0x106f8e810>]






    

In [44]:

    

plot(gus[gus.size/2-g.size:gus.size/2+1])
plot(gus2[gus2.size/2-g2.size+1:gus2.size/2+1])


    

    
Out[44]:





[<matplotlib.lines.Line2D at 0x108a7af10>]






    

In [36]:

    

plot(fft(gus[45:56]).real)
plot(gft)


    

    
Out[36]:





[<matplotlib.lines.Line2D at 0x108034f90>]






    

In [45]:

    

gft


    

    
Out[45]:





array([  1.59413217e+00 +0.00000000e+00j,
        -1.39813625e+00 -4.54282005e-01j,
         9.29589826e-01 +6.75386542e-01j,
        -4.37972931e-01 -6.02818024e-01j,
         1.09453270e-01 +3.36862526e-01j,
        -3.60822483e-16 -9.43689571e-16j,
         1.09453270e-01 -3.36862526e-01j,
        -4.37972931e-01 +6.02818024e-01j,
         9.29589826e-01 -6.75386542e-01j,  -1.39813625e+00 +4.54282005e-01j])

In [46]:

    

conj(gft)


    

    
Out[46]:





array([  1.59413217e+00 -0.00000000e+00j,
        -1.39813625e+00 +4.54282005e-01j,
         9.29589826e-01 -6.75386542e-01j,
        -4.37972931e-01 +6.02818024e-01j,
         1.09453270e-01 -3.36862526e-01j,
        -3.60822483e-16 +9.43689571e-16j,
         1.09453270e-01 +3.36862526e-01j,
        -4.37972931e-01 -6.02818024e-01j,
         9.29589826e-01 +6.75386542e-01j,  -1.39813625e+00 -4.54282005e-01j])

In [ ]: