1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
|
# This file is part of the Astrometry.net suite.
# Licensed under a 3-clause BSD style license - see LICENSE
from __future__ import print_function
import matplotlib
if __name__ == '__main__':
matplotlib.use('Agg')
import unittest
try:
import pyfits
except ImportError:
try:
from astropy.io import fits as pyfits
except ImportError:
raise ImportError("Cannot import either pyfits or astropy.io.fits")
import numpy as np
import pylab as plt
from math import pi,sqrt
from astrometry.plot.plotstuff import *
class TestPlotstuff(unittest.TestCase):
def setUp(self):
np.random.seed(42)
def test_image_wcs(self):
'''
Round-trip an image through a WCS.
'''
W,H = 100,100
I = np.random.uniform(size=(H,W))
imfn = 'test-plotstuff-1.fits'
pyfits.writeto(imfn, I, clobber=True)
wcs = anwcs_create_box(43., 11., 1., W, H)
anwcs_rotate_wcs(wcs, 32.)
wcsfn = 'test-plotstuff-1.wcs'
anwcs_write(wcs, wcsfn)
plot = Plotstuff()
plot.outformat = PLOTSTUFF_FORMAT_PNG
plot.size = (W, H)
plot.wcs_file = wcsfn
plot.color = 'black'
plot.plot('fill')
im = plot.image
plot_image_set_wcs(im, wcsfn, 0)
im.image_high = 1.
im.image_low = 0.
plot_image_set_filename(im, imfn)
plot.plot('image')
plotfn = 'test-plotstuff-1.png'
plot.write(plotfn)
I2 = plt.imread(plotfn)
print(I2.shape)
I2 = I2[:,:,0]
plt.clf()
plt.subplot(2,2,1)
plt.imshow(I, vmin=0, vmax=1)
plt.subplot(2,2,2)
plt.imshow(I2, vmin=0, vmax=1)
plt.subplot(2,2,3)
plt.imshow(I2 - I, vmin=-1, vmax=1)
plt.savefig('test-plotstuff-1b.png')
md = np.max(np.abs((I2-I).ravel()))
print('max diff', md)
self.assertTrue(md < 1./255.)
def test_subpixel_image_wcs(self):
'''
Create a well-sampled Gaussian test image, and push it through subpixel-shifted WCSes.
'''
W,H = 100,100
X,Y = np.meshgrid(np.arange(W),np.arange(H))
s = 2.
CX,CY = 32.4,42.3
I = 1./(2.*pi*s**2) * np.exp(((X-CX)**2 + (Y-CY)**2) / (-2*s**2))
I /= I.max()
print('I sum', I.sum())
print('I max', I.max())
print('first moment', np.mean(I * X) / np.mean(I), np.mean(I * Y) / np.mean(I))
imfn = 'test-plotstuff-2.fits'
pyfits.writeto(imfn, I, clobber=True)
wcs = anwcs_create_box(43., 11., 1., W, H)
anwcs_rotate_wcs(wcs, 32.)
wcsfn = 'test-plotstuff-2.wcs'
anwcs_write(wcs, wcsfn)
plot = Plotstuff()
plot.outformat = PLOTSTUFF_FORMAT_PNG
plot.size = (W, H)
plot.wcs_file = wcsfn
plot.color = 'black'
plot.plot('fill')
im = plot.image
plot_image_set_wcs(im, wcsfn, 0)
im.image_high = I.max()
im.image_low = 0.
plot_image_set_filename(im, imfn)
plot.plot('image')
plotfn = 'test-plotstuff-2.png'
plot.write(plotfn)
I2 = plt.imread(plotfn)
print(I2.shape)
I2 = I2[:,:,0]
plt.clf()
plt.subplot(2,2,1)
plt.imshow(I, vmin=0, vmax=I.max())
plt.subplot(2,2,2)
plt.imshow(I2, vmin=0, vmax=1)
plt.subplot(2,2,3)
plt.imshow(I2 - I, vmin=-1, vmax=1)
plt.savefig('test-plotstuff-2b.png')
md = np.max(np.abs((I2-I).ravel()))
print('max diff', md)
self.assertTrue(md < 1./255.)
if __name__ == '__main__':
unittest.main()
|
