source: trunk/python/edgemarker.py@ 2760

Last change on this file since 2760 was 2635, checked in by Takeshi Nakazato, 12 years ago

New Development: No

JIRA Issue: Yes CAS-2825

Ready for Test: Yes

Interface Changes: No

What Interface Changed: Please list interface changes

Test Programs: List test programs

Put in Release Notes: Yes/No

Module(s): Module Names change impacts.

Description: Describe your changes here...

Added optional plot method.


File size: 7.7 KB
Line 
1from asap.scantable import scantable
2from asap._asap import _edgemarker
3import numpy
4import math
5
6class edgemarker:
7 """
8 The edgemarker is a helper tool to calibrate OTF observation
9 without explicit OFF scans. According to a few user-specified
10 options, the class automatically detects an edge region of the
11 map and mark integrations within this region as OFF.
12
13 The edgemarker supports raster pattern as well as other generic
14 ones (e.g. lissajous, double circle). The constructor takes
15 one boolean parameter to specify whether scan pattern is raster
16 or not. This is because that edge detection algorithms for raster
17 and others are different.
18
19 Current limitation of this class is that it cannot handle some
20 complicated observed area. Typical case is that the area has
21 clear 'dent' (e.g. a composite area consisting of two diamond-
22 shaped areas that slightly overlap). In such case, the class
23 will fail to detect such feature.
24
25 Note that the class takes a copy of input data so that input
26 data will not be overwritten. Result will be provided as a
27 separate data whose contents are essentially the same as input
28 except for that some integrations are marked as OFF.
29
30 Here is typical usage:
31
32 s = scantable( 'otf.asap', average=False )
33 marker = edgemarker( israster=False )
34 marker.setdata( s )
35 marker.setoption( fraction='15%', width=0.5 )
36 marker.mark()
37
38 # get result as scantable instance
39 s2 = marker.getresult()
40
41 # save result on disk
42 marker.save( 'otfwithoff.asap', overwrite=True )
43 """
44 def __init__( self, israster=False ):
45 """
46 Constructor.
47
48 israster -- Whether scan pattern is raster or not. Set True
49 if scan pattern is raster. Default is False.
50 """
51 self.israster = israster
52 self.marker = _edgemarker( self.israster )
53 self.st = None
54
55 def setdata( self, st ):
56 """
57 Set data to be processed.
58
59 st -- Data as scantable instance.
60 """
61 self.st = st
62 self.marker._setdata( self.st, False )
63 self.marker._examine()
64
65 def setoption( self, *args, **kwargs ):
66 """
67 Set options for edge detection. Valid options depend on
68 whether scan pattern is raster or not (i.e. constructor
69 is called with israster=True or False).
70
71 === for raster (israster=True) ===
72 fraction -- Fraction of OFF integration in each raster
73 row. Either numerical value (<1.0) or string
74 is accepted. For string, its value should be
75 'auto' or format 'xx%'. For example, '10%'
76 is same as 0.1. The 'auto' option estimates
77 number of OFFs based on t_OFF = sqrt(N) t_ON.
78 Default is 0.1.
79 npts -- Number of OFF integration in each raster row.
80 Default is -1 (use fraction).
81
82 Note that number of integrations specified by the above
83 parameters will be marked as OFF from both ends. So, twice
84 of specified number/fraction will be marked as OFF. For
85 example, if you specify fraction='10%', resultant fraction
86 of OFF integrations will be 20%.
87
88 Note also that, if both fraction and npts are specified,
89 specification by npts will come before.
90
91 === for non-raster (israster=False) ===
92 fraction -- Fraction of edge area with respect to whole
93 observed area. Either numerical value (<1.0)
94 or string is accepted. For string, its value
95 should be in 'xx%' format. For example, '10%'
96 is same as 0.1. Default is 0.1.
97 width -- Pixel width for edge detection. It should be given
98 as a fraction of the median spatial separation
99 between neighboring integrations in time. Default
100 is 0.5. In the most case, default value will be fine.
101 Larger value will cause worse result. Smaller value
102 may improve result. However, if too small value is
103 set (e.g. 1.0e-5), the algorithm may not work.
104 elongated -- Set True only if observed area is extremely
105 elongated in one direction. Default is False.
106 In most cases, default value will be fine.
107 """
108 option = {}
109 if self.israster:
110 keys = [ 'fraction', 'npts' ]
111 else:
112 keys = [ 'fraction', 'width', 'elongated' ]
113 for key in keys:
114 if kwargs.has_key( key ):
115 option[key] = kwargs[key]
116
117 if len(option) > 0:
118 self.marker._setoption( option )
119
120 def mark( self ):
121 """
122 Process data. Edge region is detected according to detection
123 parameters given by setoption(). Then, integrations within
124 edge region will be marked as OFF.
125 """
126 self.marker._detect()
127 self.marker._mark()
128
129 def getresult( self ):
130 """
131 Get result as scantable instance. Returned scantable is
132 copy of input scantable except for that some data are
133 marked as OFF as a result of edge detection and marking.
134 """
135 return scantable( self.marker._get() )
136
137 def save( self, name, overwrite=False ):
138 """
139 Save result as scantable.
140
141 name -- Name of the scantable.
142 overwrite -- Overwrite existing data if True. Default is
143 False.
144 """
145 s = self.getresult()
146 s.save( name, overwrite=overwrite )
147
148 def plot( self ):
149 """
150 """
151 from matplotlib import pylab as pl
152 from asap import selector
153 from asap._asap import srctype as st
154 pl.clf()
155
156 # result as a scantable
157 s = self.getresult()
158
159 # ON scan
160 sel = selector()
161 sel.set_types( int(st.pson) )
162 s.set_selection( sel )
163 diron = numpy.array( s.get_directionval() ).transpose()
164 diron[0] = rotate( diron[0] )
165 s.set_selection()
166 sel.reset()
167
168 # OFF scan
169 sel.set_types( int(st.psoff) )
170 s.set_selection( sel )
171 diroff = numpy.array( s.get_directionval() ).transpose()
172 diroff[0] = rotate( diroff[0] )
173 s.set_selection()
174 sel.reset()
175 del s
176 del sel
177
178 # plot
179 pl.ioff()
180 ax=pl.axes()
181 ax.set_aspect(1.0)
182 pl.plot( diron[0], diron[1], '.', color='blue', label='ON' )
183 pl.plot( diroff[0], diroff[1], '.', color='green', label='OFF' )
184 [xmin,xmax,ymin,ymax] = pl.axis()
185 pl.axis([xmax,xmin,ymin,ymax])
186 pl.legend(loc='best',prop={'size':'small'},numpoints=1)
187 pl.xlabel( 'R.A. [rad]' )
188 pl.ylabel( 'Declination [rad]' )
189 pl.title( 'edgemarker result' )
190 pl.ion()
191 pl.draw()
192
193def _0to2pi( v ):
194 return v % (2.0*math.pi)
195
196def quadrant( v ):
197 vl = _0to2pi( v )
198 base = 0.5 * math.pi
199 return int( vl / base )
200
201def quadrantList( a ):
202 n = len(a)
203 nquad = numpy.zeros( 4, dtype=int )
204 for i in xrange(n):
205 v = quadrant( a[i] )
206 nquad[v] += 1
207 #print nquad
208 return nquad
209
210def rotate( v ):
211 a = numpy.zeros( len(v), dtype=float )
212 for i in xrange(len(v)):
213 a[i] = _0to2pi( v[i] )
214 nquad = quadrantList( a )
215 quadList = [[],[],[],[]]
216 rot = numpy.zeros( 4, dtype=bool )
217 if all( nquad==0 ):
218 print 'no data'
219 elif all( nquad>0 ):
220 #print 'extends in all quadrants'
221 pass
222 elif nquad[0]>0 and nquad[3]>0:
223 #print 'need rotation'
224 rot[3] = True
225 rot[2] = nquad[1]==0 and nquad[2]>0
226 #print rot
227
228 for i in xrange(len(a)):
229 if rot[quadrant(a[i])]:
230 a[i] -= 2*math.pi
231 return a
Note: See TracBrowser for help on using the repository browser.