Coverage for /wheeldirectory/casa-6.7.0-12-py3.10.el8/lib/py/lib/python3.10/site-packages/casatools/image.py: 54%
222 statements
« prev ^ index » next coverage.py v7.6.4, created at 2024-10-31 17:39 +0000
1##################### generated by xml-casa (v2) from image.xml #####################
2##################### 3ea432f56e8db3bb8bb80f239a929b71 ##############################
3from __future__ import absolute_import
4from .__casac__.image import image as _image
6from .errors import create_error_string
7from .typecheck import CasaValidator as _validator
8_pc = _validator( )
9from .coercetype import coerce as _coerce
10from .coordsys import coordsys as _wrap_coordsys
11from .componentlist import componentlist as _wrap_componentlist
12_wrap_image = lambda swig_object: image(swig_object=swig_object)
14class image:
15 _info_group_ = """images"""
16 _info_desc_ = """Operations on images"""
17 ### self
18 def __init__(self, *args, **kwargs):
19 """
20 """
21 self._swigobj = kwargs.get('swig_object',None)
22 if self._swigobj is None:
23 self._swigobj = _image()
25 def newimage(self, infile):
26 """This method is identical to ia.newimagefromfile(). The description of how it works is in the help for that method.
27 """
28 return _wrap_image(swig_object=self._swigobj.newimage(infile))
30 def newimagefromfile(self, infile):
31 """This method returns an image analysis tool associated with the specified image.
32 Constructing a image analysis tool in addition to the default ia tool allows the
33 user to operate on multiple images simultaneously. All ia.newimagefrom*()
34 methods share this behavior.
36 The parameter infile may refer to a CASA image, a Miriad image, or a FITS image.
37 FITS images of types Float are supported.
39 When finished with the newly created tool, the user should close it to free
40 up system resources (eg memory).
42 ia.newimage() is an alias for this method.
43 """
44 return _wrap_image(swig_object=self._swigobj.newimagefromfile(infile))
46 def imagecalc(self, outfile='', pixels='', overwrite=False, imagemd='', prec='float'):
47 """This method is used to evaluate a mathematical expression involving
48 existing images. It fully supports float, double, and complex float, and complex
49 double valued images.
50 The syntax of the expression supplied via the pixels
51 parameter (in what is called the Lattice Expression Language, or LEL) is
52 explained in detail in htmladdnormallink{note
53 223}{http://aips2.nrao.edu/docs/notes/223/223.html}. This is a rich
54 mathematical language with allows all manner of mathematical operations
55 to be applied to images.
57 Any image files embedded in the expression may be native casa or
58 fits (but not yet Miriad) images.
60 If successful, this method always returns an image analysis tool that
61 references the image resulting from the calculation. This returned tool
62 should always be captured and closed as soon as the user is done with it
63 to free up system resources (eg, memory). The image analysis tool on which
64 the method is called (eg the ia tool when one runs ia.imagecalc()) remains
65 unaltered, eg it still refers to the same image it did prior to the imagecalc()
66 call.
68 Values of the returned tool are evaluated "on demand". That is, only when a method
69 is run on the returned tool are the necessary values computed. And in fact, the
70 values have to be reevaluated for each operation (method call). This means that
71 there is a small performance hit for using the returned tool rather than the image
72 written to disk and that none of the images which were used in the expression
73 should be deleted while the returned tool is in use because they must be accessed
74 for calculating the expression each time an operation of the returned tool is
75 performed. These limitations do not apply to the ouput image if one is specified
76 with the outfile parameter; it is a genuine CASA image with
77 numerical values. If outfile is blank, no ouput image is written (although
78 the resulting image can still be accessed via the returned image analysis tool
79 as described below).
81 Normally you should just write the image, close the returned
82 tool, and open the results image with the default ia tool and operate on it. If
83 you are interested in conserving disk space, you don't need to keep the result of
84 the calculation around for very long, and/or you are only going to do a small
85 number of operations on the result image, should you set outfile="".
87 Note that when multiple image are used in the expression, there is
88 no garauntee about which of those images will be used to create the metadata
89 of the output image, unless imagemd is specified. If imagemd is specified, the following
90 rules of metadata copying will be followed:
92 1. The pixel data type of the image specified by imagemd and the output image must
93 be the same.
94 2. The metadata copied include the coordinate system (and so of course the dimensionality of
95 the output image must correspond to the coordinate system to be copied), the image_info record
96 (which contains things like the beam(s)), the misc_info record (should one exist in the image
97 specified by imagemd), and the units.
98 3. If the output image is a spectral image, the brightness units are set to the empty string.
99 4. If the ouptut image is a polarization angle image, the brightness unit is set to "deg" and
100 the stokes coordinate is set to have a single plane of type of Pangle.
102 The precision (float or double) of the output image pixels is determined by
103 the value specified by the prec parameter. The domain (real or complex) of
104 the output pixel values is determined from the expression.
107 """
108 return _wrap_image(swig_object=self._swigobj.imagecalc(outfile, pixels, overwrite, imagemd, prec))
110 def collapse(self, function='', axes=[ ], outfile='', region='', box='', chans='', stokes='', mask='', overwrite=False, stretch=False):
111 """This method collapses an image along a specified axis or set of axes of length N pixels to a single pixel on each
112 specified axis. Both float valued and complex valued images are supported. It computes a user-specified aggregate
113 function for pixel values along the specified axes, and places those values in the single remaining plane of
114 those axes in the output image. The method returns an image analysis tool containing the newly-created collapsed
115 image. Valid choices of aggregate functions are: 'flux' (see below for constraitns), 'madm' (median absolute
116 deviation from the median), 'max', 'mean', 'median', 'min', 'npts', 'rms', 'stddev', 'sum', 'variance' and 'xmadm'
117 (median absolute deviation from the median multipied by x, where x is the reciprocal of Phi^-1(3/4), where Phi^-1
118 is the reciprocal of the quantile function. Numerically, x = 1.482602218505602. See, eg,
119 https://en.wikipedia.org/wiki/Median_absolute_deviation#Relation_to_standard_deviation). Minimal unique matching is
120 supported for the function parameter (e.g. function = 'r' will compute the rms of the pixel values, 'med' will
121 compute the median, etc.).
123 If one specifies function='flux', the following constraints must be true:
125 1. The image must have a direction coordinate,
126 2. The image must have at least one beam,
127 3. The specified axes must be exactly the direction coordinate axes,
128 4. Only one of the non-directional axes may be non-degenerate,
129 5. The iamge brightness unit must be conformant with x*yJy/beam, where x is an optional unit (such as km/s for moments images)
130 and y is an optional SI prefix.
132 Axes may be specified as a single integer or an array of integers indicating the zero-based
133 axes along which to collapse the image. Axes may also be specified as a single or array of strings which
134 minimally and uniquely match (ignoring case) world axis names in the image (e.g. 'dec' for
135 collapsing along the declination axis or ['ri', 'd'] for collapsing along both the right ascension and
136 declination axes).
138 If outfile is not specified (or contains only whitespace characters), no image is written but the
139 collapsed image is still accessible via the image analysis tool this method always returns (which
140 references the collapsed image). If the returned object is not wanted, it should still be
141 captured and destroyed via its done() method. If this is not done, there is no guarantee
142 as to when the Python garbage collector will delete it. If the returned object is wanted, it
143 should still be deleted as soon as possible for the same reasons, e.g.
145 collapsed_image = ia.collapse(...)
146 begin{verbatim}
147 # do things (or not) with the collapsed_image and when finished working with the object, do
148 end{verbatim}
149 collapsed_image.done()
151 The reference pixel of the collapsed axis is set to 0 and its reference value is set to the mean
152 of the the first and last values of that axis in the specified region of the input image. The
153 reference value is the world coordinate value of the reference pixel. For instance, if an axis
154 to be collapsed were to be the frequency axis, in the collapsed image, the reference value would
155 be the mean value of the frequency range spanned, and would be stored in pixel 0.
157 If the input image has per plane beams, the beam at the origin of the subimage determined by
158 the selected region is arbitrarily made the global beam of the output image. In general, the user
159 should understand the pitfalls of collapsing images with multiple beams (i.e. that employing an
160 aggregate function on pixels with varying beam sizes more often than not leads to ill-defined
161 results). Convolution to a common beam is not performed automatically as part of the preprocessing
162 before the actual rebinning occurs. In such cases, therefore, the user should probably first
163 convolve the input image with a common restoring beam so that each plane has the same resolution,
164 and/or use imsmooth to smooth the data to have the same beam.
167 """
168 return _wrap_image(swig_object=self._swigobj.collapse(function, axes, outfile, region, box, chans, stokes, mask, overwrite, stretch))
170 def decimate(self, outfile='', axis=int(0), factor=int(1), method='copy', region='', mask='', overwrite=False, stretch=False):
171 """This application removes planes along the specified axis of an image. It supports both float valued and complex
172 valued images. The factor parameter represents the factor by which to reduce the number
173 of planes.
175 The method parameter represents how to calculate the pixel values of the output image. A
176 value of method="copy" means that every factorth plane of the selected region in the input
177 image will be directly copied to the corresponding plane in the output image. So, if one
178 wanted to copy every third spectral plane in the input image to the output image, one would
179 specify factor=3 and method="copy". If the selected region along the specified axis had 11
180 planes, then there would be 4 output planes which would map to planes 0, 3, 6, and 9 of
181 the specified region of input image. A value of method="mean" indicates that each of
182 factor number of planes in the range starting at each factorth plane should be averaged to
183 produce the corresponding output plane. So, if one specified factor=3 and method="mean" along
184 an axis of the selected region of the input image which had 11 pixels, the corresponding axis
185 in the output image would have three pixels and the pixel values for each of those output
186 planes would corresponding to averaging along that axis planes 0-2, 3-5, and 6-8 of the
187 selected region of the input image. Note that the remaining planes, 9 and 10, in the selected
188 region of the input image would be ignored because the last interval must have exactly
189 factor number of planes in order to be included in the output image.
191 The coordinate system of the output image takes into account the decimation; that is, along the
192 decimated axis, the increment of the output image is factor times that of the input image, and
193 the reference pixel of the output image is located at pixel 1/factor times the reference pixel
194 in the input image.
196 This method returns an image analysis tool which references the output image. If this tool
197 is not desired, one should capture it anyway and then close() it immediately to free up
198 resources.
200 Images with multiple beams are not supported; please convolve a multi-beam image to a single
201 resolution before running this application.
203 """
204 return _wrap_image(swig_object=self._swigobj.decimate(outfile, axis, factor, method, region, mask, overwrite, stretch))
206 def dohistory(self, enable=True):
207 """This allows control over if tool methods record history of what parameters they were called
208 with to the input and/or output image. By default, tool methods will write to the image
209 history. By explicitly disabling history writing, tool methods will not write to the history.
210 When created, an ia tool will have history writing enabled. Note that the setting is specific
211 to the individual tool, and that methods such as open(), close(), done(), fromshape(), etc do
212 not implicitly change the internal state of whether or not history writing by methods is
213 enabled. One can explicitly enable/disable history writing even if the tool is not yet attached
214 to an image. In the case where a method returns a new ia tool, the method will not write history
215 to the output image, but the returned tool will have history writing enabled, so that running
216 methods on the returned tool will cause history to be written to the image attached to that
217 tool, or any new image created by running methods on that tool, unless dohistory(False) is
218 explicitly run on that tool prior to running other methods.
220 IMPORTANT NOTE: This setting will not affect the behavior of ia.sethistory(); this tool
221 method will always write to the history, no matter if ia.dohistory(False) was run prior.
223 """
224 return self._swigobj.dohistory(enable)
226 def imageconcat(self, outfile='', infiles=[ ], axis=int(-1), relax=False, tempclose=True, overwrite=False, reorder=False, mode='paged'):
227 """This application is used to concatenate two or more input CASA
228 images into one output image. For example, if you have two image cubes
229 which are contiguous along one axis (say a spectral axis) and you would
230 like to glue them together along this axis, then this function is the
231 appropriate thing to use.
233 The axis parameter is used to specify which zero-based axis the images
234 should be concatenated along. A negative value indicates that the
235 spectral axis should be used. If a negative value is given but there is no
236 spectral axis, an exception will be thrown. The zero-based order of the
237 axes of an image can be determined from ia.coordsys().names().
239 If successful, this method will return an image analysis tool referencing
240 the concatenated image. Even if it is not wanted, the returned tool should
241 be captured and closed as soon as the user is finished with it to free up
242 system resources (eg memory).
244 If outfile is given, the image is written to the specified
245 disk file. If outfile is unset, the on-the-fly image tool
246 created by the method actually references all of the input files.
247 So if you deleted any of the input image disk files, it would render
248 this tool useless. Note that a blank outfile name may only be given
249 in the case where mode='p' (see below). When you destroy this tool
250 done() method, the reference connections are broken.
252 The input and output images must be of the same dimensionality. Therefore,
253 if you wish to concatenate 2-D images into a 3-D image, the 2-D images
254 must have a third axis (of length unity) so that the output image
255 coordinates are known along the concatenation axis.
257 The input images are concatenated in the order in which they are listed unless
258 the reorder parameter is set to True. If True, the images are reordered if necessary
259 so that the world coordinate values along the selected axis monotonically increase
260 or decrease. The direction of the increment is determined by the first listed image.
261 If reorder=True, the world coordinate ranges of the images along the selected axis
262 are not permitted to overlap, and the signs of the increments for this axis in all
263 images must be the same. If reorder=False, the coordinate system of the first listed
264 image is used as the coordinate system for the output image. If reorder=True, the
265 coordinate system of the first image in the list of the reordered images is used
266 as the coordinate system of the output image. Setting reorder=True can be especially
267 useful if the infiles are specified using a wildcard character(s).
269 If relax=False, the input images are checked to see that they are
270 contiguous along the concatenation axis and an exception is thrown if
271 they are not. In addition, the coordinate descriptors (e.g. reference
272 pixel, reference value etc) for the non-concatenation axes must be the
273 same or an error will result.
275 The input disk image files may be in native (CASA or FITS) formats.
277 The contiguous criterion and coordinate descriptor equality criteria can
278 be relaxed by setting relax=True whereupon only warnings will be
279 issued. Dimension and shape must still be the same though. When the
280 concatenation axis is not contiguous (but still monotonically increasing
281 or decreasing) and relax=True, a tabular coordinate will be used
282 to correctly describe the axis. But be aware that it means adjacent
283 pixels are not regularly spaced. However, methods like toworld() and
284 topixel() will correctly interconvert world and pixel coordinates.
286 In giving the input image names, the infiles argument can be a
287 single string if you wild card it with standard shell symbols. For
288 example, infiles='cena_???.*', where the '?' represents one
289 character and '*' any number of characters.
291 Otherwise, you must input a vector of strings such as
292 infiles="cena1 cena2 cena3". An input such as infiles='files1,file2'
293 will be interpreted as one string naming one file and an exception will
294 be thrown. The reason for this is that although the latter could be
295 parsed to extract two file names by recognizing comma delimiters, it
296 is not possible because an expression such as infiles='cena.{a,b}'
297 (meaning files of name 'cena.a' and 'cena.b') would confuse such
298 parsing (you would get two files of name cena.{a} and {sff b}.
300 You can look at the coordinate system of the output image using the
301 ia.summary() method to ensure it's correct.
303 The parameter tempclose is, by default, True. This means that
304 all internal reference copies of the input images are kept closed until
305 they are needed. Then they are opened temporarily and then closed again.
306 This enables you to effectively concatenate as many images as you like
307 without encountering any operating system open file number limits.
308 However, it comes at some performance loss, because opening and closing
309 all those files takes time. If you are concatenating a smallish number
310 of files, you might use tempclose=False. This will leave all internal
311 reference copies permanently open, but performance, if you don't hit the
312 file limit, will be better.
314 The mode parameter controls the format of the output image. If mode='p',
315 then a "paged" (normal persistent CASA format) image in which all the
316 pixel values and metadata are copied to a single image is created. All
317 the other modes will create a "virtual concat" image. In this format,
318 a directory, named outfile, is created to store metadata about each image.
319 The metadata are placed in a file in this directory named imageconcat.json,
320 which is a plain text file in JSON format. In the json file are pointers
321 to the actual disk images that make up the virtual concat image. In the
322 case of mode='c' (copy), alll the disk images are copied to the outfile
323 directory and the json file references the images in the outfile directory.
324 In the case of mode='m' (move), the disk images are instead moved to the
325 outfile directory, and again the json file references the images in the
326 outfile directory. These two modes make it simple to deal with the
327 virtual concat image as a single unit, for example, if the user needs
328 to move the image or tar it up, the resulting moved image or untarred file
329 will be accessible by the image tool as a valid image. The final option,
330 mode='n' (no copy) means the output directory only contains the json file
331 and references the input image paths as given. In this case, moving the
332 outfile directory without moving the reference images will normally cause
333 problems, especially if the input images were given using relative path
334 names. So, one must think carefully if mode='n' will be the correct option
335 for their use case in the long run.
337 This method requires multiple images which are specified with the infiles
338 parameter. Therefore calling ia.open() is not necessary, although calling
339 imageconcat() using an already open image analysis tool will work and the
340 state of that tool (eg the image it references) will not be changed.
341 """
342 return _wrap_image(swig_object=self._swigobj.imageconcat(outfile, infiles, axis, relax, tempclose, overwrite, reorder, mode))
344 def fromarray(self, outfile='', pixels=[ ], csys={ }, linear=False, overwrite=False, log=True, type='f'):
345 """This application converts a numerical numpy array of any size and dimensionality
346 into a CASA image. It will create float, double, complex-float, and complex-double
347 valued images.
349 The image analysis tool on which this method is called will reference the created
350 image; if this tool referenced another image before this call, that image will no
351 longer be referenced by the tool after the creation of the new image. If you
352 would rather have a new image analysis tool returned, keeping the one on which
353 this method is called unaltered, use newimagefromarray() instead. If outfile is
354 specified, a persistent image is written to disk, if not, the image tool on
355 which this method was called will reference a temporary image (either in memory
356 or on disk, depending on its size) that will be deleted when the tool is closed.
358 The type parameter controls the data type/precision of the pixel values of the
359 created image. 'f' indicates that float precision point (32 bit precision) pixel
360 values should be writted. 'd' indicates that double precision (64 bit precision)
361 pixel values should be written. If the input array has complex (as opposed to
362 real) values, then complex pixel values, with each of the real and imaginary
363 parts having the specified precision, will be written. Array values will be cast
364 automatically to the specified precision, so that the precision of the input
365 array values may be increased, decreased, or unchanged depending on the input
366 array type.
368 The coordinate system, provided as a a dictionary (use eg, cs.torecord() to do
369 that), is optional. If specified, it must have the same number of dimensions
370 as the pixels array. Call the naxes() method on the coordinate system tool to
371 see how many dimensions the coordinate system has. A coordinate system can be
372 created from scratch using the coordinate system (cs) tool and methods therein,
373 but often users prefer to use a coordinate system from an already existing image.
374 This can be gotten using ia.coordsys() which returns a coordinate system tool.
375 A torecord() call on that tool will result in a python dictionary describing
376 the coordinate system which is the necessary format for the csys input parameter
377 of ia.fromarray().
379 If csys is not specified, a default coordinate system will be created. If
380 linear=False (the default), the created coordinate system will have standard
381 RA/DEC/Stokes/Spectral Coordinate axes depending upon the shape of the pixels
382 array (Stokes axis must be no longer than 4 pixels and the spectral axis may
383 precede the Stokes axis if eg, shape=[64,64,32,4]. Extra dimensions are given
384 linear coordinates. If linear=True, then all the resulting coordinates
385 are linear with the axes represent lengths. In this case each axis will have a
386 value of 0.0 at its center pixel. The increment of each axis will be 1.0 km.
388 The method returns True if creation of the image was successful, False otherwise,
389 so a check can be made programmatically if the image creation was successful.
390 """
391 return self._swigobj.fromarray(outfile, pixels, csys, linear, overwrite, log, type)
393 def fromcomplist(self, outfile='', shape=[ ], cl='', csys={ }, overwrite=False, log=True, cache=True):
394 """This method allows one to create an image based on a component list. A component list
395 is a list of simple models (point sources, Gaussians, disks, etc) that describe the
396 sky brightness (cf the component list (cl) tool). Images that can be described in this
397 way normally require significantly less space to store than traditional images in which
398 all the pixel values are stored. For a component list image, pixel values are computed
399 "on the fly". Pixel values can be cached by specifying cache=True (the default value)
400 while the image is attached to an image tool, which permits faster access to them after
401 they are computed the first time. The trade off to caching is that resources such as
402 memory and disk space must be used to cache the pixel values.
404 The image is constrained to have two, three, or four dimensions. One must specify an
405 image shape (the dimensionality of which must adhere to this constraint). One may also
406 supply a coordinate system specification using the csys parameter. If a coordinate system
407 is not specified, a default coordinate system is used. If specified, the coordinate system
408 must have a direction coordinate which has two pixel axes. It can also have a spectral
409 and/or polarization coordinate. The maximum length of the polarization coordinate is four
410 pixels, and the world coordinate values of the polarization coordinate are constrained to
411 be in the set of stokes parameters I, Q, U, and V.
413 As is common with image creation methods, specifying an empty string for the outfile parameter
414 results in a tempoary image being created that will be deleted when either the done() or
415 close() method is run on the tool. By specifying a non-empty string, the image is saved to
416 disk and can be opened with the open() method later. A persistent component list image is
417 composed of a component list table that has metadata describing the image-related information,
418 such as the coordinate system and the shape, as well as a history (log) table.
420 Because pixel values are computed from the component models, altering pixel values is not
421 supported. So methods such as putchunk(), putregion(), and addnoise() will fail on component
422 list images when trying to modify pixel values. However, persistent image masks and on the fly
423 masks are fully supported.
425 The brightness unit of a component list image is constrained to be "Jy/pixel". Attempts to
426 modify this value using setbrightnessunit() will fail.
428 Component list images do not support synthesized beams; attempting to run setrestoringbeam() on
429 a component list image to add a beam(s) will fail.
431 One can easily create an image in which the pixel values are persistently stored from a
432 component list image by running methods such as fromimage(), subimage(), tofits(), etc. In
433 general, any method run on a component list image that creates a new image will create a
434 non-component list image (eg, a traditional CASA Paged image or Temporary image) in which the
435 pixel values are explicitly stored.
437 DISK MODELS
439 Pixels with centers inside the disk will have the same values, even if a pixel straddles the
440 edge of the disk. Pixels with straddle the edge of the disk which have centers outside the
441 disk are given values of zero. Thus, one should not expect the flux density of the disk to
442 be exactly the provided value to the component list; for a given size disk, the computed flux
443 density will be closer to the expected value for images with smaller pixels.
446 """
447 return self._swigobj.fromcomplist(outfile, shape, cl, csys, overwrite, log, cache)
449 def fromfits(self, outfile='', infile='', whichrep=int(0), whichhdu=int(0), zeroblanks=False, overwrite=False):
450 """This function is used to convert a FITS disk image file (Float,
451 Double, Short, Long are supported) to an
452 casa imagefile. If {stfaf outfile} is given, the image is written
453 to the specified disk file. If {stfaf outfile} is unset, the Image
454 tool is associated with a temporary image. This temporary image may
455 be in memory or on disk, depending on its size. When you close the
456 Image tool (with the close function) this
457 temporary image is deleted.
459 This function reads from the FITS primary array (when the image is at
460 the beginning of the FITS file; {stfaf whichhdu=0}), or an image
461 extension (when the image is elsewhere in the FITS file, {stfaf
462 whichhdu $>$ 0}).
464 By default, any blanked pixels will be converted to a mask value which
465 is false, and a pixel value that is NaN. If you set {stfaf
466 zeroblanks=T} then the pixel value will be zero rather than NaN. The
467 mask will still be set to false. See the function
468 replacemaskedpixels if you
469 need to replace masked pixel values after you have created the image.
470 """
471 return self._swigobj.fromfits(outfile, infile, whichrep, whichhdu, zeroblanks, overwrite)
473 def fromimage(self, outfile='', infile='', region=[ ], mask='', dropdeg=False, overwrite=False):
474 """This function applies a region to an imagefile, creates a new
475 imagefile containing the (sub)image, and associates the imagetool
476 with it.
478 The input image file may be in native casa, fits, or Miriad
479 format. Look htmlref{here}{IMAGES:FOREIGNIMAGES} for more
480 information on foreign images.
482 If {stfaf outfile} is given, the (sub)image is written to the specified
483 disk file.
485 If {stfaf outfile} is unset, the Image tool actually references
486 the input image file. So if you deleted the input image disk file, it
487 would render this tool useless. When you close this tool
488 (with the close function)
489 the reference connection is broken.
491 Sometimes it is useful to drop axes of length one (degenerate axes).
492 Use the {stfaf dropdeg} argument if you want to do this.
494 The output mask is the combination (logical OR) of the default input
495 pixelmask (if any) and the OTF mask. Any other input pixelmasks
496 will not be copied. Use function
497 maskhandler if you need to copy other
498 masks too.
500 See also the subimage function.
501 """
502 return self._swigobj.fromimage(outfile, infile, region, mask, dropdeg, overwrite)
504 def fitsheader(self, retstr=False, exclude=[ ]):
505 """This method constructs a FITS header dictionary which describes the current image. The
506 exclude parameter can either be a string or a list of strings. It indicates FITS header
507 keywords to exclude. This may be useful for the HISTORY keyword because the HISTORY can
508 be very long. If the retstr parameter is set to True, then a string will be returned
509 representing the FITS header. The exculde parameter has no effect when the header
510 information is returned as a single string. An new IMAGENME keyword is added to the
511 header. The value for IMAGENME is the string returned by name method, but the name is
512 truncated to the maximum size for a FITS header value, 68 characters. The image name
513 is truncated from the beginning of the string to reduce it to 68 characters, so the
514 value will contain either the full name or the final 68 characters of the name.
515 """
516 return self._swigobj.fitsheader(retstr, exclude)
518 def fromshape(self, outfile='', shape=[ int(0) ], csys={ }, linear=False, overwrite=False, log=True, type='f'):
519 """This function creates a CASA image with the specified shape. All the pixel
520 values in the image are set to 0. One may create an image with float precision
521 pixels (type='f'), complex float precision pixels (type='c'), double precision
522 pixels (type='d'), or complex double precision pixels ('cd'). To use a numpy
523 array of values to create an image, use ia.fromarray(). To make a 2-D image from
524 a packaged FITS file, use ia.maketestimage().
526 If outfile is given, the image is written to the specified disk file. If
527 outfile is unset, the image analysis tool is associated with a temporary image.
528 This temporary image may be in memory or on disk, depending on its size. When
529 you close the image analysis tool (with the ia.close() method, the temporary
530 image is deleted.
532 The coordinate system, provided as a coordinate system tool record, is optional.
533 If provided, it must be dimensionally consistent with the specified shape.
535 If the coordinate system is not provided, a default coordinate system will be
536 created. If linear=False (the default), then it is a
537 standard RA/DEC/Stokes/Spectral coordinate system depending exactly upon the
538 shape (the Stokes axis must be no longer than 4 pixels and spectral axis may
539 occur prior to the Stokes axis if eg, shape=[64,64,32,4]. Extra dimensions are
540 given linear coordinates. If linear=True, then the coordinate system will have
541 linear coordinates.
543 The method returns True if successful, False otherwise.
544 """
545 return self._swigobj.fromshape(outfile, shape, csys, linear, overwrite, log, type)
547 def maketestimage(self, outfile='', overwrite=False):
548 """This function converts a FITS file resident in the casa system into
549 a casa image.
551 If outfile is given, the image is written to the specified disk
552 file. If outfile is unset, the Image tool is associated with a
553 temporary image. This temporary image may be in memory or on disk,
554 depending on its size. When you close the Image tool (with the close
555 function) this temporary image is deleted.
556 """
557 return self._swigobj.maketestimage(outfile, overwrite)
559 def deviation(self, outfile='', region='', mask='', overwrite=False, stretch=False, grid=[ int(1),int(1) ], anchor='ref', xlength='1pix', ylength='1pix', interp='cubic', stattype='sigma', statalg='classic', zscore=float(-1), maxiter=int(-1)):
560 """This application creates an image that reflects the statistics of the input image. The output image has
561 the same dimensions and coordinate system as the (selected region in) input image. The grid parameter
562 describes how many pixels apart "grid" pixels are. Statistics are computed around each grid pixel. Grid
563 pixels are limited to the direction plane only; independent statistics are computed for each direction plane
564 (ie at each frequency/stokes pixel should the input image happen to have such additional axes). Using the
565 xlength and ylength parameters, one may specify either a rectangular or circular region around each grid
566 point that defines which surrounding pixels are used in the statistic computation for individual grid points.
567 If the ylength parameter is the empty string, then a circle of diameter provided by xlength centered on
568 the grid point is used. If ylength is not empty, then a rectangular box of dimensions xlength x ylength centered
569 on the grid pixel is used. These two parameters may be specified in pixels, using either numerical values or
570 valid quantities with "pix" as the unit (eg "4pix"). Otherwise, they must be specified as valid angular
571 quantities, with recognized units (eg "4arcsec"). As with other region selections in CASA, full pixels are
572 included in the computation even if the specified region includes only a fraction of that pixel. BEWARE OF
573 MACHINE PRECISION ISSUES, because you may get a smaller number of pixels included in a region than you
574 expect if you specify, eg, an integer number of pixels. In such cases, you probably want to specify that
575 number plus a small epsilon value (eg "2.0001pix" rather than "2pix") to mitigate machine precision issues
576 when computing region extents.
578 The output image is formed by putting the statistics calculated at each grid point at the corresponding
579 grid point in the output image. Interpolation of these output values is then used to compute values at
580 non-grid-point pixels. The user may specify which interpolation algorithm to use for this computation
581 using the interp parameter.
583 The input image pixel mask is copied to the output image. If interpolation is performed, output pixels are
584 masked where the interpolation fails.
586 ANCHORING THE GRID
588 The user may choose at which pixel to "anchor" the grid. For example, if one specifies grid=[4,4] and
589 anchor=[0,0], grid points will be located at pixels [0,0], [0,4], [0,8] ... [4,0], [4,4], etc. This
590 is exactly the same grid that would be produced if the user specified anchor=[4,4] or anchor=[20,44].
591 If the user specifies anchor=[1, 2] and grid=[4,4], then the grid points will be at pixels [1,2], [5,2],
592 [9,2]... [5,2], [5,6], etc. and the resulting grid is the same as it would be if the user specified eg
593 anchor=[9,10] or anchor=[21, 18]. The value "ref", which is the default, indicates that the reference
594 pixel of the input image should be used to anchor the grid. The x and y values of this pixel will be
595 rounded to the nearest integer if necessary.
597 SUPPORTED STATISTICS AND STATISTICS ALGORITHMS
599 One may specify which statistic should be represented using the stattype parameter. The following values
600 are recognized (minimum match supported):
602 +--------------------+-------------------------------------------------------------------------------------------------------+
603 | iqr | inner quartile range (q3 - q1) |
604 +--------------------+-------------------------------------------------------------------------------------------------------+
605 | max | maximum |
606 +--------------------+-------------------------------------------------------------------------------------------------------+
607 | mean | mean |
608 +--------------------+-------------------------------------------------------------------------------------------------------+
609 | medabsdevmed, madm | median absolute deviation from the median |
610 +--------------------+-------------------------------------------------------------------------------------------------------+
611 | median | median |
612 +--------------------+-------------------------------------------------------------------------------------------------------+
613 | min | minimum |
614 +--------------------+-------------------------------------------------------------------------------------------------------+
615 | npts | number of points |
616 +--------------------+-------------------------------------------------------------------------------------------------------+
617 | q1 | first quartile |
618 +--------------------+-------------------------------------------------------------------------------------------------------+
619 | q3 | third quartile |
620 +--------------------+-------------------------------------------------------------------------------------------------------+
621 | rms | rms |
622 +--------------------+-------------------------------------------------------------------------------------------------------+
623 | sigma, std | standard deviation |
624 +--------------------+-------------------------------------------------------------------------------------------------------+
625 | sumsq | sum of squares |
626 +--------------------+-------------------------------------------------------------------------------------------------------+
627 | sum | sum |
628 +--------------------+-------------------------------------------------------------------------------------------------------+
629 | var | variance |
630 +--------------------+-------------------------------------------------------------------------------------------------------+
631 | xmadm | median absolute deviation from the median multipied by x, where x is the reciprocal of Phi^-1(3/4), |
632 | | where Phi^-1 is the reciprocal of the quantile function. Numerically, x = 1.482602218505602. See, eg, |
633 | | https://en.wikipedia.org/wiki/Median_absolute_deviation#Relation_to_standard_deviation |
634 +--------------------+-------------------------------------------------------------------------------------------------------+
636 Using the statalg parameter, one may also select whether to use the Classical or Chauvenet/ZScore statistics algorithm to
637 compute the desired statistic (see the help for ia.statistics() or imstat for a full description of these algorithms).
640 """
641 return _wrap_image(swig_object=self._swigobj.deviation(outfile, region, mask, overwrite, stretch, grid, anchor, xlength, ylength, interp, stattype, statalg, zscore, maxiter))
643 def adddegaxes(self, outfile='', direction=False, spectral=False, stokes='', linear=False, tabular=False, overwrite=False, silent=False):
644 """This method adds degenerate axes (i.e.
645 axes of length 1) of the specified type. Sometimes this can be useful
646 although you will generally need to modify the coordinate system of the
647 added axis to give it the coordinate you want (do this with the
648 Coordsys tool). This method supports
649 both float and complex valued images.
651 You specify which type of axes you want to add. You can't add
652 an axis type that already exists in the image. For the Stokes axis,
653 the allowed value (a string such as I, Q, XX, RR) can be found in the
654 Coordsys newcoordsys function documentation.
656 If {stfaf outfile} is given, the image is written to the specified
657 disk file. If {stfaf outfile} is unset, the on-the-fly Image tool
658 returned by the function is associated with a temporary image. This
659 temporary image may be in memory or on disk, depending on its size.
660 When you destroy the generated Image tool (with the done function) this
661 temporary image is deleted.
662 """
663 return _wrap_image(swig_object=self._swigobj.adddegaxes(outfile, direction, spectral, stokes, linear, tabular, overwrite, silent))
665 def addnoise(self, type='normal', pars=[ float(0.0),float(1.0) ], region={ }, zero=False, seeds=[ ]):
666 """This function adds noise to the image. You may zero the image first
667 before the noise is added if you wish.
669 The noise can be drawn from one of many distributions.
671 For each distribution, you must supply the type via the {stfaf type}
672 argument (minimum match is active) and parameters via the {stfaf
673 pars} argument. Briefly:
675 begin{itemize}
677 item {binomial} -- the binomial distribution models successfully drawing
678 items from a pool. Specify two parameters, $n$ and $p$, respectively.
679 $n$ is the number of items in the pool, and $p$, is the probability of
680 each item being successfully drawn. It is required that $n > 0$ and
681 $0 le p le 1$.
683 item {discreteuniform} -- models a uniform random variable over the closed interval. Specify
684 two parameters, the low and high values, respectively.
685 The low parameter is the lowest possible return value and
686 the high parameter is the highest. It is required that $low < high$.
688 item {erlang} -- Specify two parameters, the mean and variance,
689 respectively. It is required that the mean is non-zero and the variance
690 is positive.
692 item {geometric} -- Specify one parameter, the probability.
693 It is required that $0 le probability < 1$.
695 item {hypergeometric} -- Specify two parameters, the mean and the variance.
696 It is required that the variance is positive and that the mean is non-zero
697 and not bigger than the square-root of the variance.
699 item {normal} -- Specify two parameters, the mean and the variance.
700 It is required that the variance is positive.
702 item {lognormal} -- Specify two parameters, the mean and the variance.
703 It is required that the supplied variance is positive and that the mean is non-zero.
705 item {negativeexponential} -- Supply one parameter, the mean.
707 item {poisson} -- Specify one parameter, the mean.
708 It is required that the mean is non-negative.
710 item {uniform} -- Model a uniform random variable over a closed
711 interval. Specify two parameters, the low and high values. The low
712 parameter is the lowest possible return value and the high parameter can
713 never be returned. It is required that $low < high$.
715 item {weibull} -- Specify two parameters, alpha and beta.
716 It is required that the alpha parameter is not zero.
718 The random number generator seeds may be specified as an array of integers. Only the first
719 two values are used. If none or a single value is provided, the necessary remaining value(s)
720 are generated based on the current time, using the algorithm
721 begin{verbatim}
722 seedBase = 1e7*MJD
723 seed[1] = (Int)seedBase;
724 # and if seed[0] is also not supplied
725 seed[0] = (Int)((1e7*(seedBase - seed[1])))
726 end{verbatim}
728 where MJD is the Modidfied Julian Day.
730 end{itemize}
731 """
732 return self._swigobj.addnoise(type, pars, region, zero, seeds)
734 def convolve(self, outfile='', kernel=[ ], scale=float(-1.0), region={ }, mask='', overwrite=False, stretch=False):
735 """This function performs Fourier-based convolution of an image by the specified
736 kernel. The input image must have real-valued pixels.
738 If outfile is specified, a persistent image is written to the specified disk
739 file. If outfile is left set to the empty string, the on-the-fly image analysis
740 tool generated by this function is associated with a temporary image. This
741 temporary image may be stored in memory or on disk, depending on its size.
742 When the user destroys the generated image analysis tool with the close() or
743 done() method, the temporary image is deleted.
745 The kernel is provided as a multi-dimensional array or as the filename of a
746 persistent image. The provided kernel can have fewer dimensions than the image
747 being convolved. In this case, it will be padded with degenerate axes. An
748 exception will be thrown if the kernel has more dimensions than the image. No
749 additional scaling of the kernel is provided.
751 The scaling of the output image is determined by the value of the scale
752 parameter. If this is left unset, then the kernel is normalized to unit sum.
753 If scale is not left unset, then the convolution kernel will be scaled
754 (multiplied) by this value.
756 Masked pixels will be assigned the value 0.0 before convolution.
758 The output mask is the logical AND of the input image's default pixel mask (if
759 any) and the OTF mask. Any other pixel masks associated with the input image
760 will not be copied. The maskhandler method should be used if there is a need to
761 copy other masks too.
763 See also the other convolution methods convolve2d(), sepconvolve(), and
764 hanning() for more specialized convolution algorithms.
765 """
766 return _wrap_image(swig_object=self._swigobj.convolve(outfile, kernel, scale, region, mask, overwrite, stretch))
768 def boundingbox(self, region={ }):
769 """This function finds the bounding box of a
770 region of interest when it is applied to a particular image. Both
771 float and complex valued images are supported. It is
772 returned in a record which has fields {cf 'blc', 'trc', 'inc',
773 'bbShape', 'regionShape', 'imageShape', 'blcf'} and {cf 'trcf'}
774 containing the bottom-left corner, the top-right corner (in absolute
775 image pixel coordinates), the increment (stride) of the region, the
776 shape of the boundingbox, the shape of the region, the shape of the
777 image, the blc in formatted absolute world coordinates and the trc in
778 formatted absolute world coordinates, respectively.
780 Note that the shape of the bounding box will be different from the shape
781 of the region if a non-unit stride (increment) is involved (see the example
782 below).
784 Note that the integer size of the elements in blc, trc, inc, regionShape,
785 bbShape, and imageShape are 32 bits, even on a 64 bit machine. This means that,
786 on 64 bit machines, you may have to convert them to 64 bit ints using eg
787 numpy.int64, before being able to use them as direct input to other
788 methods such as ia.getchunk().
789 """
790 return self._swigobj.boundingbox(region)
792 def boxcar(self, outfile='', region={ }, mask=[ ], axis=int(-1), width=int(2), drop=True, dmethod='copy', overwrite=False, stretch=False):
793 """This application performs boxcar convolution of one axis of an image
794 defined by
796 z[i] = (y[i] + y[i+i] + ... + y[i+w-1])/w
798 where z[i] is the value at pixel i in the box car smoothed image, y[k]
799 is the pixel value of the input image at pixel k, and w is a postivie integer
800 representing the width of the boxcar in pixels. Both float and complex
801 valued images are supported. The length of the axis along which the
802 convolution is to occur must be at least w pixels in the selected region,
803 unless decimation using the mean function is chosen in which case the axis
804 length must be at least 2*w (see below). Masked pixel values are set to
805 zero prior to convolution. All nondefault pixel masks are ignored during
806 the calculation. The convolution is done in the image domain (i.e., not
807 with an FFT).
809 If drop=False (no decimation), the length of the output axis will be equal
810 to the length of the input axis - w + 1. The pixel mask, ORed with the OTF mask
811 if specified, is copied from the selected region of the input image to the
812 output image. Thus for example, if the selected region in the input image has
813 six planes along the convolution axis, if the specified boxcar width is 2,
814 and if the pixel values, which are all unmasked, on a slice along this axis
815 are [1, 2, 5, 10, 17, 26], then the corresponding output slice will be of
816 length five and the output pixel values will be [1.5, 3.5, 7.5, 13.5, 21.5].
818 If drop=True and dmethod="copy", the output image is the image calculated
819 if drop=True, except that only every wth plane is kept. Both the pixel and mask
820 values of these planes are copied directly to the output image, without further
821 processing. Thus for example, if the selected region in the input image has six
822 planes along the convolution axis, the boxcar width is chosen to be 2, and if
823 the pixel values, which are all unmasked, on a slice along this axis are [1, 2,
824 5, 10, 17, 26], the corresponding output pixel values will be [1.5, 7.5, 21.5].
826 If drop=True and dmethod="mean", first the image described in the drop=False
827 case is calculated. Then, the ith plane of the output image is calculated by
828 averaging the i*w to the (i+1)*w-1 planes of this intermediate image. Thus, for
829 example, if the selected region in the input image has six planes along the
830 convolution axis, the boxcar width is chosen to be 2, and if the pixel values,
831 which are all unmasked, on a slice along this axis are [1, 2, 5, 10, 17, 26],
832 then the corresponding output pixel values will be [2.5, 10.5]. Any pixels at the
833 end of the plane of the intermediate image that do not fall into a complete bin of
834 width w are ignored. Masked values are taken into consideration when forming this
835 average, so if one of the values is masked, it is not used in the average. If at
836 least one of the valuesin the intermediate image bin is not masked, the
837 corresponding output pixel will not be masked.
839 The smoothed image is written to disk with name {stfaf outfile}, if specified.
840 If not, no image is written but the image is still accessible via the returned
841 image analysis tool (see below).
843 This method always returns an image analysis tool which is attached to the smoothed
844 image. This tool should always be captured and closed after any desired manipulations
845 have been done. Closing the tool frees up system resources (eg memory), eg,
847 begin{verbatim}
848 smoothedim = ia.boxcar(...)
849 # do things (or not) with smoothedim
850 ...
851 # close the returned tool promptly upon finishing with it.
852 smoothedim.done()
853 end{verbatim}
855 """
856 return _wrap_image(swig_object=self._swigobj.boxcar(outfile, region, mask, axis, width, drop, dmethod, overwrite, stretch))
858 def brightnessunit(self):
859 """This function gets the image brightness unit. Both float and complex
860 valued images are supported.
861 """
862 return self._swigobj.brightnessunit()
864 def calc(self, pixels, verbose=True):
865 """This function is used to evaluate a mathematical expression involving
866 casa images, assigning the result to the current (already existing)
867 image. Both float and complex valued images are supported, although the
868 image which results from the calculation must have the same type of pixel
869 values as the image already attached to the tool. That is, one cannot
870 create a complex valued image using this method if the associated ia tool
871 is currently attached to a float valued image. It complements the imagecalc
872 function which returns a newly constructed on-the-fly image tool. See htmladdnormallink{note 223}{../../notes/223/223.html}
873 which describes the the syntax and functionality in detail.
875 If the expression, supplied via the {stfaf pixels} argument, is not a
876 scalar, the shapes and coordinates of the image and expression must
877 conform.
879 If the image (that associated with the tool) has a pixelmask, then only
880 pixels for which that mask is good will be changed. See the function
881 maskhandler for managing image pixelmasks.
883 Note that when multiple image are used in the expression, there is no garauntee about which of
884 those images will be used to create the header of the output image. Therefore, one may have
885 to modify the output header as needed if the input headers differ.
887 See the related functions set and
888 putregion.
889 """
890 return self._swigobj.calc(pixels, verbose)
892 def calcmask(self, mask='', name='', asdefault=True):
893 """This method is used to create a new pixel mask via a boolean LEL expression.
895 See http://casa.nrao.edu/aips2_docs/notes/223/index.shtml which describes the
896 the syntax and functionality of LEL in detail. Also in this document is a
897 description of ways to escape image names that contain certain non-alphanumeric
898 characters so they are compatible with LEL syntax.
900 If the expression is not a scalar, the shapes and coordinates of the image and
901 expression must conform. If the expression is a scalar then the entire pixel
902 mask will be set to that value.
904 By specifying the name parameter to be an empty string, the method automatically
905 determines the name of the output mask. Otherwise, the output mask is named the
906 value specified by the name parameter. If a mask by that name already exists,
907 it is overwritten. Use method ia.summary() to view existing mask names.
909 The asdefault parameter specifies if the new mask should be set as the default
910 pixel mask of the image. This value is set to True by default.
911 """
912 return self._swigobj.calcmask(mask, name, asdefault)
914 def close(self):
915 """This function closes the imagetool. This means that it detaches the
916 tool from its imagefile (flushing all the changes first). The
917 imagetool is ``null'' after this change (it is not destroyed) and
918 calling any toolfunction other than open
919 will result in an error.
920 """
921 return self._swigobj.close()
923 def continuumsub(self, outline='', outcont='continuumsub.im', region={ }, channels=[ int(-1) ], pol='', fitorder=int(0), overwrite=False):
924 """This function packages the relevant image tool functionality for simple
925 specification and application of image plane continuum subtraction. All
926 that is required of the input image is that it have a non-degenerate
927 spectral axis.
929 The user specifies region, the region of the input image over which
930 continuum subtraction is desired (otherwise the whole image will be
931 treated); channels, the subset of channels on the spectral axis to use
932 in the continuum estimation, specified as a vector;
933 fitorder, the polynomial order to use in the
934 estimation. Optionally, output line and continuum images may be written
935 by specifying outline and outcont, respectively. If outline is not
936 specified, a virtual image tool is all that is produced. If outcont is
937 not specified, the output continuum image will be written in
938 'continuumsub.im'. Note that the pol parameter is no longer supported; one
939 should use the region parameter if polarization selection is desired, in
940 conformance with other ia tool methods.
941 """
942 return _wrap_image(swig_object=self._swigobj.continuumsub(outline, outcont, region, channels, pol, fitorder, overwrite))
944 def convertflux(self, value=[ ], major=[ ], minor=[ ], type='Gaussian', topeak=True, channel=int(-1), polarization=int(-1)):
945 """This function interconverts between peak intensity and flux density for a
946 Gaussian component. The image must hold a restoring beam.
947 """
948 return self._swigobj.convertflux(value, major, minor, type, topeak, channel, polarization)
950 def convolve2d(self, outfile='', axes=[ int(0),int(1) ], type='gaussian', major='0deg', minor='0deg', pa='0deg', scale=float(-1), region={ }, mask='', overwrite=False, stretch=False, targetres=False, beam={ }):
951 """This function performs Fourier-based convolution of an imagefile
952 using the provided 2D kernel.
954 If {stfaf outfile} is left unset, the image is written to the specified
955 disk file. If {stfaf outfile} is not given, the newly constructed
956 on-the-fly Image tool is associated with a temporary image. This
957 temporary image may be stored in memory or on disk, depending on its size.
958 When the user destroys the on-the-fly Image tool (with the done function) this
959 temporary image is deleted.
961 The user specifies which 2 pixel axes of the image are to be convolved
962 via the {stfaf axes} argument. The pixels must be square or an error will result.
964 The user specifies the type of convolution kernel with {stfaf type} (minimum
965 match is supported); currently only {cf 'gaussian'} is available.
967 The user specifies the parameters of the convolution kernel via the arguments
968 {stfaf major}, {stfaf minor}, and {stfaf pa}. These arguments can
969 be specified in one of three ways:
971 begin{itemize}
973 item Quantity - for example {stfaf major=qa.quantity(1, 'arcsec')}
974 Note that you pixel units can be used, viz. {stfaf major=qa.quantity(1, 'pix')},
975 see below.
977 item String - for example {stfaf minor='1km'} (i.e. one that the
978 Quanta quantity function accepts).
980 item Numeric - for example {stfaf major=10}. In this case, the units
981 of {stfaf major} and {stfaf minor} are assumed to be in pixels. Using
982 pixel units allows the user to convolve unlike axes (see one of the provided
983 example for this use case). For the position angle, units of degrees are assumed.
985 end{itemize}
987 The interpretation of {stfaf major} and {stfaf minor} depends upon the
988 kernel type.
990 begin{itemize}
992 item Gaussian - {stfaf major} and {stfaf minor} are
993 the Full Width at Half Maximum (FWHM) of the major and minor
994 axes of the Gaussian.
996 end{itemize}
998 The position angle is measured North through East when a
999 plane holding a celestial coordinate (the usual astronomical
1000 convention) is convolved. For other axis/coordinate combinations,
1001 a positive position angle is measured from +x to +y in the
1002 absolute pixel coordinate frame (x is the first axis that is
1003 specified, with argument {stfaf axes}).
1005 In the case of a Gaussian, the {stfaf beam} parameter offers an alternate way of
1006 describing the convolving Gaussian. If used, neither {stfaf major}, {stfaf minor},
1007 nor {stfaf pa} can be specified. The {stfaf beam} parameter must have exactly three
1008 fields: "major", "minor", and "pa" (or "positionangle"). This is, not coincidentally,
1009 the record format for the output of ia.restoringbeam().
1011 The scaling of the output image is determined by the argument {stfaf scale}.
1012 If this is left unset then autoscaling will be invoked.
1014 If the user is not convolving the sky, then autoscaling means that the convolution
1015 kernel will be normalized to have unit volume so as to conserve flux.
1017 If the user is convolving the sky, then there are two cases
1018 for which autoscaling is useful:
1020 Firstly, if the input image units are Jy/pixel, then the output image
1021 will have units of Jy/beam and be appropriately scaled. In addition,
1022 the restoring beam of the output image will be the same as the
1023 convolution kernel.
1025 Secondly,if the input image units are Jy/beam, then the output
1026 image will also have units of Jy/beam and be appropriately
1027 scaled. In addition, the restoring beam of the output image
1028 will be the convolution of the input image restoring beam and the
1029 convolution kernel. In the case of an image with per-plane beams, for
1030 each plane, the kernel is convolved with the appropriate beam and the
1031 result is associated with that plane in the output image.
1033 If the user sets a value for {stfaf scale}, then the convolution kernel
1034 will be scaled by this value. Note that it has peak of unity before the
1035 application of this scale factor.
1037 If the units on the original image include Jy/beam, the units on the
1038 output image will be rescaled by the ratio of the input and output
1039 beams as well as rescaling by the area of convolution kernel.
1041 If the units on the original image include K, then only the image
1042 convolution kernel rescaling is done.
1044 If targetres=True and type="gaussian" and the input image has a restoring beam,
1045 this method will interpret the values of major, minor, and pa as the resolution
1046 of the final image and will calculate the parameters of the Gaussian to use
1047 in the convolution so that this target resolution is achieved.
1049 Masked pixels will be assigned the value 0.0 before convolution.
1050 The output mask is the intersection (logical AND) of the default input
1051 pixelmask (if any) and the OTF mask. Any other input pixelmasks
1052 will not be copied. The function
1053 maskhandler
1054 can be used if there is a need to copy other masks too.
1056 See also the other convolution functions:
1058 convolve,
1059 hanning, and
1060 sepconvolve.
1061 """
1062 return _wrap_image(swig_object=self._swigobj.convolve2d(outfile, axes, type, major, minor, pa, scale, region, mask, overwrite, stretch, targetres, beam))
1064 def coordsys(self, axes=[ int(-1) ]):
1065 """This function returns the Coordinate System of an image in a {stf
1066 Coordsys} tool. Both float and complex valued images are supported.
1068 By default, the Coordinate System describes all of the axes in the
1069 image. If you desire, you can select a subset of the axes, thus
1070 reducing the dimensionality of the Coordinate System. This may be
1071 useful if you are supplying a Coordinate System to the
1072 functions fromarray or
1073 fromshape.
1074 """
1075 return _wrap_coordsys(swig_object=self._swigobj.coordsys(axes))
1077 def coordmeasures(self, pixel=[ float(-1) ], dframe='cl', sframe='cl'):
1078 """You can use this function to get the world coordinates for a specified
1079 absolute pixel coordinate in the image. You specify a pixel coordinate
1080 (0-rel) for each axis in the image.
1082 If you supply fewer pixel values then there are axes in the image, your
1083 value will be padded out with the reference pixel for the missing axes.
1084 Excess values will be ignored.
1086 The parameters dframe and sframe allow one to specify to which reference frame
1087 the direction and spectral measures, respectively, should be converted. These
1088 values are case-insensitive. "native" means use the native reference frame of
1089 the coordinate in question. "cl" means use the conversion layer frame if one
1090 exists (if not, the native frame will be used).
1092 The world coordinate is returned as a record of measures. This
1093 function is just a wrapper for the Coordsys tool toworld function
1094 (invoked with argument {stfaf format='m'}). Please see its
1095 documentation for discussion about the formatting and meaning of the
1096 measures.
1098 This Image tool function adds two additional fields to the return record.
1100 The {cf mask} field contains the value of the image pixelmask at the
1101 specified position. It is either T (pixel is good) or F (pixel is masked
1102 as bad or the specified position was off the image).
1104 The {cf intensity} field contains the value of the image (at the
1105 nearest pixel to that given) and its units. This is actually stored
1106 as a Quantity. This field does not exist
1107 if the specified pixel coordinate is off the image.
1108 """
1109 return self._swigobj.coordmeasures(pixel, dframe, sframe)
1111 def decompose(self, region={ }, mask='', simple=False, threshold=float(-1), ncontour=int(11), minrange=int(1), naxis=int(2), fit=True, maxrms=float(-1), maxretry=int(-1), maxiter=int(256), convcriteria=float(0.0001), stretch=False):
1112 """This function is an image decomposition tool that performs several tasks,
1113 with the end result being that a strongly blended image is separated into
1114 components - both in the sense that it determines the parameters for each
1115 component (assuming a Gaussian model) and that it physically assigns each
1116 pixel in the image to an individual object. The products of these two
1117 operations are called the component list and the component map,
1118 respectively. The fitting process (which determines the component list) and
1119 the pixel-decomposition process (which determines the component map) are
1120 designed to work cooperatively to increase the efficiency and accuracy of
1121 both.
1123 The algorithm behind the decomposition is based on the function clfind,
1124 described in Williams et al 1994, which uses a contouring procedure whereby
1125 a closed contour designates a separate component. The program first
1126 separates the image into clearly distint 'regions' of blended emission, then
1127 contours each region to determine the areas constituting each component and
1128 passes this information on to the fitter, which determines the component
1129 list.
1131 The contour deblending can optionally be replaced with a simpler local maximum
1132 scan, and the fitting can be replaced with a moment-based estimation method to
1133 speed up calculations on very large images or if either primary method causes
1134 trouble, but in general this will impede the accuracy of the fit.
1136 The function works with both two and three dimensional images.
1138 The return value is a record (or dictionary) that has 3 keys: {tt 'components', 'blc', 'trc'}.
1139 The {tt 'components'} element is a matrix each row of which contains the gaussian parameters of the component fitted.
1140 The {tt 'blc'} element is a matrix of the bottom left corners (blc) of the regions found. Each row correspond to a region blc.
1141 The {tt 'trc'} element is a matrix of the top right corners (trc) of the regions found. Each row correspond to a region trc.
1142 {bf Please Note} that the returned blc's and trc's are relative to {tt region} defined by the user. A {tt blc } of [0,0] implies the bottom left of the region selected and not the bottom left of the image. Obviously if no region is defined then it is the bottom left of the image.
1144 """
1145 return self._swigobj.decompose(region, mask, simple, threshold, ncontour, minrange, naxis, fit, maxrms, maxretry, maxiter, convcriteria, stretch)
1147 def deconvolvecomponentlist(self, complist, channel=int(-1), polarization=int(-1)):
1148 """This method deconvolves (a
1149 record representation of) a Componentlist tool from the restoring
1150 beam, returning (a record representation of) a new Componentlist tool.
1151 If there is no restoring beam, a fail is generated.
1153 Currently, only deconvolution of Gaussian components is supported.
1155 For images with per-plane beam, the user must choose which beam is used for
1156 the deconvolution by setting channel and/or polarization. Only a single beam
1157 is used to deconvolve all components.
1159 See also functions setrestoringbeam and
1160 restoringbeam.
1161 """
1162 return self._swigobj.deconvolvecomponentlist(complist, channel, polarization)
1164 def deconvolvefrombeam(self, source=[ ], beam=[ ]):
1165 """This is a helper function. It is to provide a way to deconvolve gaussians from other gaussians if that is what is needed for example removing a beam Gaussian from a Gaussian source. To run this function the tool need not be attached to an image.
1167 The return value is a record that contains the fit param and the return value is a boolean which is set to true if fit model is a point source
1168 """
1169 return self._swigobj.deconvolvefrombeam(source, beam)
1171 def beamforconvolvedsize(self, source=[ ], convolved=[ ]):
1172 """Determine the size of the beam necessary to convolve with the given source to reach the
1173 given convolved (source+beam) size. Because the problem is completely specified by the
1174 input parameters, no image needs to be attached to the associated tool; eg ia.open() need
1175 not be called prior to calling this method.
1176 """
1177 return self._swigobj.beamforconvolvedsize(source, convolved)
1179 def commonbeam(self):
1180 """Determine a beam to which all beams in an image can be convolved.
1181 If the image does not have a beam, an exception will be thrown.
1182 If the image has a single beam, that beam will be returned.
1183 If the image has multiple beams, this will be the beam with the largest area in the image
1184 beam set if all the other beams can be convolved to that beam. If not, this is guaranteed to be the minimum area beam to which
1185 all beams in the set can be convolved if all but one of the beams in the set can be convolved to the beam in the set with the
1186 largest area. Otherwise, the returned beam may or may not be the smallest possible beam to which all the beams in the set
1187 can be convolved.
1189 """
1190 return self._swigobj.commonbeam()
1192 def remove(self, done=False, verbose=True):
1193 """This function first closes the
1194 imagetool which detaches it from its underlying imagefile. It then
1195 deletes that imagefile. If {stfaf done=False}, the imagetool is still
1196 viable, and can be used with function open
1197 to open a new imagefile. Otherwise the imagetool is destroyed. If {stfaf verbose=True}, the logger will receive a progress report.
1198 """
1199 return self._swigobj.remove(done, verbose)
1201 def removefile(self, file):
1202 """This function deletes the specified image file.
1203 """
1204 return self._swigobj.removefile(file)
1206 def done(self, remove=False, verbose=True):
1207 """When the user no longer needs to use an imagetool, calling this function
1208 will free up its resources. That is, it destroys the tool. This means
1209 that the user can no longer call any functions on the tool after it
1210 has been {stff done}.
1212 If the Image tool is associated with a disk file, then (unlike the
1213 {stff close} function, the user can also choose to delete that by
1214 setting {stfaf remove=true}. By default, any associated disk file is
1215 not deleted.
1217 Note that this function is different from the {stff close} function
1218 because the latter does not destroy the imagetool. For example, the
1219 user can use the {stff open} function straight after the {stff close}
1220 function on the same tool.
1221 """
1222 return self._swigobj.done(remove, verbose)
1224 def fft(self, real='', imag='', amp='', phase='', axes=[ int(-1) ], region={ }, mask='', stretch=False, complex=''):
1225 """This method fast Fourier Transforms the supplied image to the Fourier plane.
1226 If the axes parameter is left unset, then the direction plane of the image (if
1227 there is one) is transformed. Otherwise, the user can specify which axes are
1228 to be transformed. Note that if the direction plane is to be transformed, both
1229 axes associated with it must be specified.
1231 The user specifies which form is desired in the result by specifying the
1232 desired output image file name(s).
1234 Before the FFT is performed, any masked pixels are set to values of zero. The
1235 output mask is the result of ANDing the default input pixel mask (if any) and
1236 the OTF mask. Any other input pixel masks will not be copied. The method
1237 maskhandler() can be used if there is a need to copy other masks.
1239 The following rules are used to set the brightness units of the output images.
1240 1. The phase image always has units of radians.
1241 For the other output images,
1242 2. if the input image has units of Jy/beam or Jy/pixel
1243 (ie it is an image-plane image), the output (uv-plane) images will have units
1244 of Jy. In the case of the input image having a synthesized beam, the beam
1245 will be copied to the output images (which is important for transforming back).
1246 3. If the input image has units of Jy (ie is a uv-plane image), the output
1247 images will have either units of Jy/beam or Jy/pixel, depending on if the
1248 input image has a beam or not.
1250 For some transformations (e.g., UV domain to image domain transforms), it is
1251 not possible to automatically generate an expected coordinate system for
1252 the output image(s); only the FFT numerics are performed and the coordinate
1253 system is generated using generic conventions.
1254 """
1255 return self._swigobj.fft(real, imag, amp, phase, axes, region, mask, stretch, complex)
1257 def findsources(self, nmax=int(20), cutoff=float(0.1), region={ }, mask='', point=True, width=int(5), negfind=False):
1258 """This function finds strong point sources in
1259 the image. The sources are returned in a record that can be used by a
1260 Componentlist tool.
1262 An efficient method is used to locate sources under the assumption that
1263 they are point-like and not too close to the noise. Only sources with a
1264 peak greater than the {stfaf cutoff} fraction of the strongest source
1265 will be found. Only positive sources will be found, unless the {stfaf
1266 negfind=T} whereupon positive and negative sources will be found.
1268 After the list of point sources has been made, you may choose to make a
1269 Gaussian fit for each one ({stfaf point=F}) so that shape information
1270 can be recovered as well. You can specify the half-width of the
1271 fitting grid with argument {stfaf width} which defaults to 5 (fitting
1272 grid would then be [11,11] pixels). If you set {stfaf width=0}, this is
1273 a signal that you would still like Gaussian components returned, but a
1274 default width should be used for the Gaussian shapes. The default is
1275 such that the component is circular with a FWHM of {stfaf width}
1276 pixels.
1278 Thus, if {stfaf point=T}, the components in the returned Componentlist
1279 are Point components. If {stfaf point=F} then Gaussian components are
1280 returned.
1282 The region must be 2-dimensional and it must hold a region of the sky.
1283 Any degenerate trailing dimensions in the region are discarded.
1285 See also the function fitcomponents (for which {stff
1286 findsources} can provide an initial estimate).
1287 """
1288 return self._swigobj.findsources(nmax, cutoff, region, mask, point, width, negfind)
1290 def fitprofile(self, box='', region='', chans='', stokes='', axis=int(-1), mask='', ngauss=int(1), poly=int(-1), estimates='', minpts=int(1), multifit=False, model='', residual='', amp='', amperr='', center='', centererr='', fwhm='', fwhmerr='', integral='', integralerr='', stretch=False, logresults=True, pampest=[ ], pcenterest=[ ], pfwhmest=[ ], pfix='', gmncomps=int(0), gmampcon=[ ], gmcentercon=[ ], gmfwhmcon=[ ], gmampest=[ float(0.0) ], gmcenterest=[ float(0.0) ], gmfwhmest=[ float(0.0) ], gmfix='', spxtype='', spxest=[ ], spxfix=[ ], div=[ ], spxsol='', spxerr='', logfile='', append=True, pfunc='', goodamprange=[ float(0.0) ], goodcenterrange=[ float(0.0) ], goodfwhmrange=[ float(0.0) ], sigma='', outsigma='', planes=[ ]):
1291 """This application simultaneously fits any number of gaussian singlets, any number of lorentzian singlets, and any number of gaussian multiplets,
1292 and/or a polynomial to one dimensional profiles using the non-linear, least squares Levenberg-Marquardt algorithm. A description of the
1293 fitting algorithm may be found in AIPS++ Note 224 (http://www.astron.nl/casacore/trunk/casacore/doc/notes/224.html) and in Numerical Recipes
1294 by W.H. Press et al., Cambridge University Press. A gaussian/lorentzian singlet is a gaussian/lorentzian whose parameters (amplitude,
1295 center position, and width) are all independent from any other feature that may be simultaneously fit. A gaussian multiplet is a set of two or
1296 more gaussian lines in which at least one (and possibly two or three) parameter of each line is dependent on the parameter of another,
1297 single (reference) profile in the multiplet. For example, one can specify a doublet in which the amplitude of the first line is 0.6 times the
1298 amplitude of the zeroth line and/or the center of the first line is 20 pixels from the center of the zeroth line, and/or the fwhm of the first
1299 line is identical (in pixels) to that of the zeroth line. There is no limit to the number of components one can specify in a multiplet
1300 (except of course that the number of parameters to be fit should be significantly less than the number of data points), but there can be only
1301 a single reference profile in a multiplet to which to tie constraints of parameters of the other profiles in the set.
1303 Additionally, a power logarithmic polynomial (plp) or a logarithmic tranformed polynomial (ltp) can be fit. In this case, each of these functions
1304 cannot be fit simultaneously with any other supported function. These functions are most often used for fitting the spectral index and
1305 higher order terms of a spectrum. A power logarithmic polynomial has the form
1307 y = c0*x/div**(c1 + c2*ln(x/div) + c3*ln(x/div)**2 + ... + cn*ln(x/div)**(n - 1))
1309 and a logarithmic transformed polynomial is simply the result of this equation after taking the natural log of both sides so that it has the form
1311 ln(y) = c0 + c1*ln(x/div) + c2*ln(x/div)**2 + ... + cn*ln(x/div)**n
1313 The coefficients of the two forms correspond with each other except that c0 in the second equation is equal to
1314 ln(c0) of the first. In the case of fitting a spectral index, the spectral index, traditionally represented as alpha, is
1315 equal to c1.
1317 In both cases, div is a numerical value used to scale abscissa values so they are closer to unity when they are sent to the fitter. This generally
1318 improves the probability that the fit will converge. This parameter may be specified via the div parameter. A value of 0
1319 (the default) indicates that the application should determine a reasonable value for div, which is determined via
1321 div = 10**int(log10(sqrt(min(x)*max(x))))
1323 where min(x) and max(x) are the minimum and maximum abscissa values, respectively.
1325 So, for example, if S(nu) is proportional to nu**alpha and you expect alpha to be near -0.8 and the value of S(nu) is 1.5 at
1326 1e9 Hz and your image(s) have spectral units of Hz, you would specify spxest=[1.5, -0.8] and div=1e9 when fitting a plp function,
1327 or spxest=[0.405, -0.8] and div=1e9 if fitting an ltp function.
1329 More details of fitting all of these functions are described in following sections.
1331 A CAUTIONARY NOTE
1332 Note that the likelihood of getting a reliable solution increases with the number of good data points as well as the goodness
1333 of the initial estimate. It is possible that the first solution found might not be the best one, and
1334 so, if a solution is found, it is recommended that the fit be repeated using the solution of the previous fit as the
1335 initial estimatE for the new fit. This process should be repeated until the solutions from one fit to the next differ only insignificantly.
1336 The convergent solution is very likely the best solution.
1338 AXIS
1339 The axis parameter indicates on which axis profiles should be fit; a value <0 indicates the spectral axis should be used, or if one does not exist,
1340 that the zeroth axis should be used.
1342 MINIMUM NUMBER OF PIXELS
1343 The minpts parameter indicates the minimum number of unmasked pixels that must be present in order for a fit
1344 to be attempted. When multifit=True, positions with too few good points will be masked in any output images.
1346 ONE FIT OF REGION AVERAGE OR PIXEL BY PIXEL FIT
1347 The multifit parameter indicates if profiles should be fit at each pixel in the selected region (true), or if the profiles in that region should be
1348 averaged and the fit done to that average profile (false).
1350 POLYNOMIAL FITTING
1351 The order of the polynomial to fit is specified only via the poly parameter. If poly<0, no polynomial will be fit. No initial estimates of
1352 coefficients can be specified; these are determined automatically.
1354 GAUSSIAN SINGLET FITTING
1355 In the absence of an estimates file and no estimates being specified by the p*est parameters, and gmncomps=0 or is empty, the ngauss parameter
1356 indicates the maximum number of gaussian singlets that should be fit. The initial estimates of the parameters for these gaussians will be attempted
1357 automatically in this case. If it deems appropriate, the fitter will fit fewer than this number. In the case where an estimates file is supplied,
1358 ngauss is ignored (see below). ngauss is also ignored if the p*est parameters are specified or if gmncomps is not an empty array or, if an integer,
1359 is greater than zero. If estimates is not specified or the p*est parameters are not specified and ngauss=0, gmncomps is empty or 0, and poly<0,
1360 an error will occur as this indicates there is nothing to fit.
1362 One can specify initial estimates of gaussian singlet parameters via an estimates file or the pampest, pcenterest, pfwhmest, and optionally, the
1363 pfix parameters. The latter is the recommended way to specify these estimates as support for estimates files may be deprecated in the future. No matter
1364 which option is used, an amplitude initial estimate must always be nonzero. A negative fwhm estimate will be silently changed to positve.
1366 SPECIFYING INITIAL ESTIMATES FOR GAUSSIAN AND LORENTZIAN SINGLETS (RECOMMENDED METHOD)
1367 One may specify initial estimates via the pampest, pcenterest, and pfwhmest parameters. In the case of a single gaussian or lorentzian singlet,
1368 these parameters can be numbers. pampest must be specified in image brightness units, pcenterest must be given in the number of pixels from the
1369 zeroth pixel, and pfwhmest must be given in pixels. Optionally pfix can be specified and in the case of a single gaussian or lorentzian singlet
1370 can be a string. In it is coded which parameters should be held constant during the fix. Any combination of "p" (amplitude), "c" (center), or "f"
1371 (fwhm) is allowed; eg pfix="pc" means fix both the amplitude and center during the fit. In the case of more than one gaussian and/or lorentzian
1372 singlets, these parameters must be specified as arrays of numbers. The length of the arrays indicates the number of singlets to fit and must be
1373 the same for all the p*est parameters.
1375 If no parameters are to be fixed for any of the singlets, pfix can be set to the empty string. However, if at least one parameter of one singlet
1376 is to be fixed, pfix must be an array of strings and have a length equal to the p*est arrays. Singlets which are not to have any parameters fixed
1377 should be represented as an empty string in the pfix array. So, for example, if one desires to fit three singlets and fix the fwhm of the middle
1378 one, one must specify pfix=["", "f", ""], the empty strings indicating no parameters of the zeroth and second singlet should be held constant.
1380 In the case of multifit=True, the initial estimates, whether from the p*est parameters or from a file (see below), will be applied to the location
1381 of the first fit. This is normally the bottom left corner of the region selected. If masked, not enough good points to perform a fit, or the
1382 attempted fit fails, the fitting proceeds to the next pixel with the pixel value of the lowest numbered axis changing the fastest. Once a
1383 successful fit has been performed, subsequent fits will use the results of a fit for a nearest pixel for which a previous fit was successful as the
1384 initial estimate for the parameters at the current location. The fixed parameter string will be honored for every fit performed when multifit=True.
1386 One specifies what type of PCF profile to fit via the pfunc parameter. A PCF function is one that can be parameterized by a peak, center, and FWHM,
1387 as both gaussian and lorentzian singlets can. If all singlets to be fit are gaussians, one can set pfunc equal to the empty string and all snglets
1388 will be assumed to be gaussians. If at least one lorentzian is to be fit, pfunc must be specified as a string (in the case of a single singlet) or
1389 an array of strings (in the case of multiple singlets). The position of each string corresponds to the positions of the initial estimates in the
1390 p*est and pfix arrays. Minimal match ("g", "G", "l", or "L") is supported. So, if one wanted to simultaneously fit two gaussian and two lorentzian
1391 singlets, the zeroth and last of which were lorentzians, one would specify pfunc=["L", "G", "G", "L"].
1393 ESTIMATES FILE FOR GAUSSIAN SINGLETS (NONRECOMMENDED METHOD)
1394 Initial estimates for gaussian singlets can be specified in an estimates file. Estimates files may be deprecated in the future in favor of the
1395 p*est parameters, so it is recommended users use those parameters instead. If an estimates file is desired to be used, the p*est parameters
1396 must be 0 or empty and mgncomps must be 0 or empty. Only gaussian singlets can be specified in an estimates file. If one desires to fit one or
1397 more gaussian multiplets and/or one or more lorentzian singlets simultaneously, the p*est parameters must be used to specify the initial parameters
1398 of all gaussian singlets to fit; one cannot use an estimates file in this case. If an estimates file is specified, a polynomial
1399 can be fit simultaneously by specifying the poly parameter. The estimates file must contain initial estimates of parameters
1400 for all gaussian singlets to be fit. The number of gaussian singlets to fit is gotten from the number of estimates in the file. The file can contain
1401 comments which are indicated by a "#" at the beginning of a line. All non-comment lines will be interpreted as initial estimates. The
1402 format of such a line is
1404 [peak intensity], [center], [fwhm], [optional fixed parameter string]
1406 The first three values are required and must be numerical values. The peak intensity must be expressed in image brightness units, while the
1407 center must be specified in pixels offset from the zeroth pixel, and fwhm must be specified in pixels. The fourth value is optional and if present,
1408 represents the parameter(s) that should be held constant during the fit. Any combination of the characters 'p' (peak), 'c' (center), and 'f' (fwhm) are
1409 permitted, eg "fc" means hold the fwhm and the center constant during the fit. Fixed parameters will have no error associated with them. Here is an
1410 example file:
1412 begin{verbatim}
1413 # estimates file indicating that two gaussians should be fit
1414 # first guassian estimate, peak=40, center at pixel number 10.5, fwhm = 5.8 pixels, all parameters allowed to vary during
1415 # fit
1416 40, 10.5, 5.8
1417 # second gaussian, peak = 4, center at pixel number 90.2, fwhm = 7.2 pixels, hold fwhm constant
1418 4, 90.2, 7.2, f
1419 # end file
1420 end{verbatim}
1422 GAUSSIAN MULTIPLET FITTING
1423 Any number of gaussian multiplets, each containing any number of two or more components, can be simultaneously fit, optionally with a
1424 polynomial and/or any number of gaussian and/or lorentzian singlets, the only caveat being that the number of parameters to be fit should be
1425 significantly less than the number of data points. The gmncomps parameter indicates the number of multiplets to fit and the number of
1426 components in each multiplet. In the case of a single multiplet, an integer (>1) can be specified. For example, mgncomps=4 means fit a
1427 single quadruplet of gaussians. In the case of 2 or more multiplets, and array of integers (all >1) must be specified. For example,
1428 gmncomps=[2, 4, 3] means 3 seperate multiples are to be fit, the zeroth being a doublet, the first being a quadruplet, and the second
1429 being a triplet.
1431 Initial estimates of all gaussians in all multiplets are specified via the gm*est parameters which must be arrays of numbers. The order
1432 starts with the zeroth component of the zeroth multiplet to the last component of the zeroth multiplet, then the zeroth component of
1433 the first multiplet to the last compoenent of the first multiplet, etc to the zeroth component of the last multiplet to the last
1434 element of the last multiplet. The zeroth element of a multiplet is defined as the reference component of that multiplet and has the special
1435 significance that it is the profile to which all constraints of all other profiles in that multiplet are referenced (see below). So,
1436 in our example of gmncomps=[2, 4, 3], gmampest, gmcenterest, and gmfwhmest must each be nine (the total number of individual gaussian
1437 profiles summed over all multiplets) element arrays. The zeroth, second, and sixth elements represent parameters of the reference profiles
1438 in the zeroth, first, and second multiplet, respectively.
1440 The fixed relationships between the non-reference profile(s) and the reference profile of a multiplet are specified via the gmampcon,
1441 gmcentercon, and gmfwhmcon parameters. At least one, and any combination, of constraints can be specified for any non-reference
1442 component of a multiplet. The amplitude ratio of a non-reference line to that of the reference line is set in gmampcon. The ratio of
1443 the fwhm of a non-reference line to that of the reference line is set in gmfwhmcon. The offset in pixels of the center position of
1444 a non-reference line to that of the reference line is set in gmcentercon. In the case where a parameter is not constrained for any
1445 non-reference line of any multiplet, the value of the associated parameter must be 0. In the case of
1446 a single doublet, a constraint may be specified as a number or an array of a single number. For example, mgncomps=2 and gmampcon=0.65
1447 and gmcentercon=[32.4] means there is a single doublet to fit where the amplitude ratio of the first to the zeroth line is constained
1448 to be 0.65 and the center of the first line is constrained to be offset by 32.4 pixels from the center of the zeroth line. In cases
1449 of a total of three or more gaussians, the constraints parameters must be specified as arrays with lengths equal to the total number
1450 of gaussians summed over all multiplets minus the number of reference lines (one per multiplet, or just number of multiplets, since
1451 reference lines cannot be constrained by themselves). In the cases where an array must be specified but a component in that array
1452 does not have that constraint, 0 should be specified. Here's an example
1454 gmncomps=[2, 4, 3]
1455 gmampcon= [ 0 , 0.2, 0 , 0.1, 4.5, 0 ]
1456 gcentercon=[24.2, 45.6, 92.7, 0 , -22.8, -33.5]
1457 gfwhmcon=""
1459 In this case we have our previous example of one doublet, one quadruplet, and one triplet. The first component of the doublet has the constraint
1460 that its center is offset by 24.2 pixels from the zeroth (reference) component. The first component of the quadruplet is constrained to have
1461 an amplitude of 0.2 times that of the quadruplet's zeroth component and its center is constrained to be offset by 45.6 pixels from the
1462 reference component. The second component of the quadruplet is constained to have its center offset by 92.7 pixels from the associated
1463 reference component and the third component is constrained to have an amplitude of 0.1 times that of the associated reference component.
1464 The first component of the triplet is constrained to have an amplitude of 4.5 times that of its associated reference component and its center
1465 is constrained to be offset by -22.8 pixels from the reference component's center. The second component of the triplet is constrained to have
1466 its center offset by -33.5 pixels from the center of the reference component. No lines have FWHM constraints, so the empty string can be given
1467 for that parameter. Note that using 0 to indicate no constraint for line center means that one cannot specify a line centered at the same
1468 position as the reference component but having a different FWHM from the reference component. If you must specify this very unusual case,
1469 try using a very small positive (or even negative) value for the center constraint.
1471 Note that when a parameter for a line is constrained, the corresponding value for that component in the corresponding gm*est array is
1472 ignored and the value of the constrained parameter is automatically used instead. So let's say, for our example above, we had specified
1473 the following estimates:
1475 gmampest = [ 1, .2, 2, .1, .1, .5, 3, 2, 5]
1476 gmcenterest = [20, 10 , 30, 45.2, 609 , -233, 30, -859, 1]
1478 Before any fitting is done, the constraints would be taken into account and these arrays would be implicitly rewritten as:
1480 gmampest = [ 1, .2, 2, .4, .1, .2, 3, 13.5, 5 ]
1481 gmcenterest = [20, 44.2, 30, 75.6, 127.7, -233, 30, 7.2, -3.5]
1483 The value of gmfwhmest would be unchanged since there are no FWHM constraints in this example.
1485 In addition to be constrained by values of the reference component, parameters of individual components can be fixed. Fixed parameters
1486 are specified via the gmfix parameter. If no parameters are to be fixed, gmfix can be specified as the empty string or a zero element
1487 array. In the case where any parameter is to be fixed, gmfix must be specified as an array of strings with length equal to the total number of
1488 components summed over all multiplets. These strings encode which parameters to be fixed for the corresponding components. If
1489 a component is to have no parameters fixed, an empty string is used. In other cases one or more of any combination of parameters can
1490 be fixed using "p", "c", and/or "f" described above for fixing singlet parameters. There are a couople of special cases
1491 to be aware of. In the case where a non-reference component parameter is constrained and the corresponding reference component parameter is
1492 set as fixed, that parameter in the non-reference parameter will automatically be fixed even if it was specified not to be fixed in
1493 the gmfix array. This is the only way the constraint can be honored afterall. In the converse case of when a constrained parameter of a
1494 non-reference component is specified as fixed, but the corresponding parameter in the reference component is not specified to be fixed,
1495 an error will occur. Fixing an unconstrained parameter in a non-reference component is always legal as is fixing any combination of
1496 parameters in a reference component (with the above caveat that corresponding constrained parameters in non-reference components will
1497 be silently held fixed as well).
1499 The same rules that apply to singlets when multifit=True apply to multiplets.
1501 LIMITING RANGES FOR SOLUTION PARAMETERS
1502 In cases of low (or no) signal to noise spectra, it is still possible for the fit to converge, but often to a
1503 nonsensical solution. The astronomer can use her knowledge of the source to filter out obviously bogus solutions.
1504 Any solution which contains a NaN value as a value or error in any one of its parameters is automatically marked as
1505 invalid.
1507 One can also limit the ranges of solution parameters to known "good" values via the goodamprange, goodcenterrange, and goodfwhmrange
1508 parameters. Any combination can be specified and the limit constraints will be ANDed together. The ranges apply to all PCF components
1509 that might be fit; choosing ranges on a component by component basis is not supported. If specified,
1510 an array of exactly two numerical values must be given to indicate the range of acceptable solution values for
1511 that parameter. goodamprange is expressed in terms of image brightness units. goodcenterrange is expressed in terms of pixels
1512 from the zeroth pixel in the specified region. goodfwhmrange is expressed in terms of pixels (only non-negative values should be
1513 given for FWHM range endpoints). In the case of a multiple-PCF fit, if any of the corresponding solutions are outside the specified
1514 ranges, the entire solution is considered to be invalid.
1516 In addition, solutions for which the absolute value of the ratio of the amplitude error to the amplitude exceeds 100 or the
1517 ratio of the FWHM error to the FWHM exceeds 100 are automatically marked as invalid.
1519 POWER LOGARITHMIC POLYNOMIAL AND LOGARITHMIC TRANSFORMED POLYNOMIAL FITTING
1520 Fitting of a sngle power logarithmic polynomial or a single logarithmic transformed polynomial function is supported.
1521 No other functions may be fit simultaneously with either of these; if parameters relating to other functions are supplied
1522 simultaneously with parameters relating
1523 to these functions, an exception will occur. For details of the functional forms, see the introduction of this
1524 document.
1526 The set of c0 ... cn coefficients (as defined previously) can
1527 be solved for. Initial estimates for the c values should be supplied via the plpest or ltpest parameters, depending on which
1528 form is being fit. The number of values given
1529 in this array will be the number of coeffecients that are solved for. One may specify which coefficients should be held
1530 fixed during the fit in the plpfix or ltpfix array. If supplied, this array should have the same number of elements as its respective
1531 initial estimates array. A value
1532 of True means the corresponding coefficient will be held fixed during the fit. An empty array indicates that no
1533 parameters will be held fixed. This is the default.
1535 Because the logarithm of the ordinate values must be taken before fitting a logarithmic transformed polynomial,
1536 all non-positive pixel values are effectively masked for the purposes of fitting.
1538 INCLUDING STANDARD DEVIATIONS OF PIXEL VALUES
1539 If the standard deviations of the pixel values in the input image are known and they vary in the image (eg they are higher for pixels
1540 near the edge of the band), they can be included in the sigma parameter. This parameter takes either an array or an image name. The
1541 array or image must have one of three shapes: 1. the shape of the input image, 2. the same dimensions as the input image with the lengths
1542 of all axes being one except for the fit axis which must have length corresponding to its length in the input image, or 3. be one
1543 dimensional with lenght equal the the length of the fit axis in the input image. In cases 2 and 3, the array or pixels in sigma will
1544 be replicated such that the image that is ultimately used is the same shape as the input image. The values of sigma must be non-negative.
1545 It is only the relative values that are important. A value of 0 means that pixel should not be used in the fit. Other than that, if pixel
1546 A has a higher standard deviation than pixel B, then pixel A is noisier than pixel B and will receive a lower weight when the fit is done.
1547 The weight of a pixel is the usual
1549 weight = 1/(sigma*sigma)
1551 In the case of multifit=False, the sigma values at each pixel along the fit axis in the hyperplane perpendicular to the fit axis which includes
1552 that pixel are averaged and the resultant averaged standard deviation spectrum is the one used in the fit. Internally, sigma values are normalized
1553 such that the maximum value is 1. This mitigates a known overflow issue.
1555 One can write the normalized standard deviation image used in the fit but specifying its name in outsigma. This image can then be
1556 used as sigma for subsequent runs.
1559 RETURNED DICTIONARY STRUCTURE
1560 The returned dictionary has a (necessarily) complex structure. First, there are keys "xUnit" and "yUnit" whose values are
1561 the abscissa unit and the ordinate unit described by simple strings. Next there are arrays giving a broad overview of the
1562 fit quality. These arrays have the shape of the specified region collapsed along the fit axis with the axis corresponding to the fit
1563 axis having length of 1:
1565 attempted: a boolean array indicating which fits were attempted (eg if too few unmasked points, a fit will not be attempted).
1566 converged: a boolean array indicating which fits converged. False if the fit was not attempted.
1567 valid: a boolean array indicating which solutions fall within the specified valid ranges of parameter space (see
1568 . section LIMITING RANGES FOR SOLUTION PARAMETERS for details).
1569 niter: an int array indicating the number of iterations for each profile, <0 if the fit did not converge
1570 ncomps: the number of components (gaussian singlets + lorentzian singlets + gaussian multiplets + polynomial) fit for the profile,
1571 . <0 if the fit did not converge
1572 direction: a string array containing the world direction coordinate for each profile
1574 There is a "type" array having number of dimensions equal to the number of dimensions in the above arrays plus one. The shape of
1575 the first n-1 dimensions is the same as the shape of the above arrays. The length of the last dimension is equal to the number of
1576 components fit. The values of this array are strings describing the components that were fit at each possition ("POLYNOMIAL",
1577 "GAUSSIAN" in the case of gaussian singlets, "LORENTZIAN" in the case of lorentzian singlets, and ""GAUSSIAN MULTPLET").
1579 If any gaussian singlets were fit, there will be a subdictionary accessible via the "gs" key which will have subkeys "amp", "ampErr", "center",
1580 "centerErr", "fwhm", "fwhmErr, "integral", and "integralErr". Each of these arrays will have one more dimension than the overview arrays described
1581 above. The shape of the first n-1 dimensions will be the same as the shape of the arrays described above, while the final dimension will
1582 have length equal to the maximum number of gaussian singlets that were fit. Along this axis will be the
1583 corresponding fit result or associated error (depending on the array's associated key) of the fit for that singlet component number. In cases where
1584 the fit did not converge, or that particular component was excluded from the fit, a value of NAN will be present.
1586 If any lorentzian singlets were fit, their solutions will be accessible via the "ls" key. These arrays follow the same rules
1587 as the "gs" arrays described above.
1589 If any gaussian multiplets were fit, there will be subdictionaries accessible by keys "gm0", "gm1", ..., "gm{n-1}" where n is the number of gaussian
1590 muliplets that were fit. Each of these dictionaries will have the same arrays described above for gaussian singlets. The last dimension
1591 will have length equal to the number of components in that particular multiplet. Each pixel along the last axis will be the parameter solution
1592 value or error for that component number in the multiplet, eg the zeroth pixel along that axis contains
1593 the parameter solution or error for the reference component of the multiplet.
1595 The polynomial coefficient solutions and errors are not returned, although they are logged.
1597 If a power logarithmic polynomial was fit, there will be a subdictionary accessible via the "plp" key which will have
1598 subkeys "soltuion" and "error" which will each have an array value. Each of these arrays will have one more dimension than the overview arrays
1599 described above. The shape of the first n-1 dimensions will be the same as the shape of the overview arrays described above, while the
1600 final dimension will have length equal to the number of parameters that were fit. Along this axis will be the
1601 corresponding fit result or associated error (depending on the array's associated key) of the fit. In cases where
1602 the fit was not attempted or did not converge, a value of NAN will be present.
1604 OUTPUT IMAGES
1605 In addition to the returned dictionary, optionally one or more of any combination of output images can be written.
1606 The model and residual parameters indicate the names of the model and residual images to be written; blank values inidcate that these images
1607 should not be written.
1609 One can also write none, any or all of the solution and error images for gaussian singlet, lorentzian singlet, and gaussian multiplet fits
1610 via the parameters amp, amperr, center, centererr, fwhm, fwhmerr, integral, integralerr when doing multi-pixel fits. For a power logarithmic
1611 polynomial or a logarithmic transformed polynomial fit, plpsol or ltpsol and plperr or ltpsol are the names of the solution and error
1612 images to write, respectively.
1614 These images contain the arrays described for the associated parameter solutions or errors described in previous sections. Each
1615 component is written to a different image, and each image is distiguished by the component it represents by its name ending
1616 in an uderscore and the relevant component number ("_0", "_1", etc). In the case of Gaussian multiplets, the image name ends
1617 with the number of the mulitplet group followed by the number of the component in that group (eg "_3_4" represents component
1618 4 in multiplet group 3). In the case of lorentzian singlets, "_ls" is appended to the image names (but before the
1619 identifying component number), in the case of gaussian multiplets. Similarly "_gm" is included in the name of Gaussian multiplet
1620 images. Pixels for which fits were not attempted, did not converge, or converged but have values of NaN (not a number) or
1621 INF (infinity) will be masked as bad.
1623 Writing analogous images for polynomial coefficients is not supported.
1626 """
1627 return self._swigobj.fitprofile(box, region, chans, stokes, axis, mask, ngauss, poly, estimates, minpts, multifit, model, residual, amp, amperr, center, centererr, fwhm, fwhmerr, integral, integralerr, stretch, logresults, pampest, pcenterest, pfwhmest, pfix, gmncomps, gmampcon, gmcentercon, gmfwhmcon, gmampest, gmcenterest, gmfwhmest, gmfix, spxtype, spxest, spxfix, div, spxsol, spxerr, logfile, append, pfunc, goodamprange, goodcenterrange, goodfwhmrange, sigma, outsigma, planes)
1629 def fitcomponents(self, box='', region=[ ], chans=[ ], stokes='', mask='', includepix=[ float(-1) ], excludepix=[ float(-1) ], residual='', model='', estimates='', logfile='', append=True, newestimates='', complist='', overwrite=False, dooff=False, offset=float(0.0), fixoffset=False, stretch=False, rms='', noisefwhm='', summary=''):
1630 """OVERVIEW
1632 This application is used to fit one or more two dimensional gaussians to sources in an image as
1633 well as an optional zero-level offset. Fitting is limited to a single polarization
1634 but can be performed over several contiguous spectral channels.
1635 If the image has a clean beam, the report and returned dictionary will contain both the convolved
1636 and the deconvolved fit results.
1638 When dooff is False, the method returns a dictionary with keys named 'converged', 'pixelsperarcsec',
1639 'results', and 'deconvolved'. The value of 'converged' is a boolean array which indicates if the
1640 fit converged on a channel by channel basis. The value of 'pixelsperarcsec' is a two element double
1641 array with the absolute values of the direction coordinate pixel increments (longitude-like and
1642 latitude-like coordinate, respectively) in arcsec. The value of 'results' is a dictionary
1643 representing a component list reflecting the fit results. In the case of an image containing beam
1644 information, the sizes and position angles in the 'results' dictionary are those of the source(s)
1645 convolved with the restoring beam, while the same parameters in the 'deconvolved' dictionary
1646 represent the source sizes deconvolved from the beam. In the case where the image does not
1647 contain a beam, 'deconvolved' will be absent. Both the 'results' and 'deconvolved' dictionaries can
1648 be read into a component list tool (default tool is named cl) using the fromrecord() method for
1649 easier inspection using tool methods, eg
1651 cl.fromrecord(res['results'])
1653 although this only works if the flux density units are conformant with Jy.
1655 There are also values in each component subdictionary not used by cl.fromrecord() but meant to
1656 supply additional information. There is a 'peak' subdictionary for each component that provides the
1657 peak intensity of the component. It is present for both 'results' and 'deconvolved' components.
1658 There is also a 'sum' subdictionary for each component indicated the simple sum of pixel values in
1659 the the original image enclosed by the fitted ellipse. There is a 'channel' entry in the 'spectrum'
1660 subdictionary which provides the zero-based channel number in the input image for which the solution
1661 applies. In addtion, if the image has a beam(s), then there will be a 'beam' subdictionary associated
1662 with each component in both the 'results' and 'deconvolved' dictionaries. This subdictionary will
1663 have three keys: 'beamarcsec' will be a subdictionary giving the beam dimensions in arcsec,
1664 'beampixels' will have the value of the beam area expressed in pixels, and 'beamster' will have the
1665 value of the beam area epressed in steradians. Also, if the image has a beam(s), in the component level
1666 dictionaries will be an 'ispoint' entry with an associated boolean value describing if the component
1667 is consistent with a point source. Each component level dictionary will have a 'pixelcoords' entry
1668 which has the value of a two element numeric array which provides the direction pixel coordinates
1669 of the fitted position.
1671 If dooff is True, in addtion to the specified number of
1672 gaussians, a zero-level offset will also be fit. The initial estimate for this
1673 offset is specified using the offset parameter. Units are assumed to be the
1674 same as the image brightness units. The zero level offset can be held constant during
1675 the fit by specifying fixoffset=True. In the case of dooff=True, the returned
1676 dictionary contains two additional keys, 'zerooff' and 'zeroofferr', which are both
1677 dictionaries containing 'unit' and 'value' keys. The values associated with the 'value'
1678 keys are arrays containing the the fitted zero level offset value and its error, respectively,
1679 for each channel. In cases where the fit did not converge, these values are set to NaN.
1680 The value associated with 'unit' is just the image brightness unit.
1682 The region can either be specified by a box(es) or a region.
1683 Ranges of pixel values can be included or excluded from the fit. If specified using
1684 the box parameter, multiple boxes can be given using the format
1685 box="blcx1, blcy1, trcx1, trcy1, blcx2, blcy2, trcx2, trcy2, ... , blcxN, blcyN, trcxN, trcyN"
1686 where N is the number of boxes. In this case, the union of the specified boxes will be used.
1688 If specified, the residual and/or model images for successful fits will be written.
1690 If an estimates file is not specified, an attempt is made to estimate
1691 initial parameters and fit a single Gaussian. If a multiple Gaussian fit
1692 is desired, the user must specify initial estimates via a text file
1693 (see below for details).
1695 The user has the option of writing the result of the fit to a log file,
1696 and has the option of either appending to or overwriting an existing file.
1698 The user has the option of writing the (convolved) parameters of a successful
1699 fit to a file which can be fed back to fitcomponents() as the estimates file for a
1700 subsequent run.
1702 The user has the option of writing the fit results in tabular format to a file whose
1703 name is specified using the summary parameter.
1705 If specified and positive, the value of rms is used to calculate the parameter uncertainties,
1706 otherwise, the rms in the selected region in the relevant channel is used for these calculations.
1708 The noisefwhm parameter represents the noise-correlation beam FWHM. If specified as a quantity,
1709 it should have angular units. If specified as a numerical value, it is set equal to that number
1710 of pixels. If specified and greater than or equal to the pixel size, it is used to calculate
1711 parameter uncertainties using the correlated noise equations (see below). If it is specified but
1712 less than a pixel width, the the uncorrelated noise equations (see below) are used to
1713 compute the parameter uncertainties. If it is not specified and the image has a restoring beam(s),
1714 the the correlated noise equations are used to compute parameter uncertainties using the
1715 geometric mean of the relevant beam major and minor axes as the noise-correlation beam FWHM. If
1716 noisefwhm is not specified and the image does not have a restoring beam, then the uncorrelated
1717 noise equations are used to compute the parameter uncertainties.
1719 SUPPORTED UNITS
1721 Currently only images with brightness units conformant with Jy/beam, Jy.km/s/beam, and K are fully
1722 supported for fitting. If your image has some other base brightness unit, that unit will be assumed
1723 to be equivalent to Jy/pixel and results will be calculated accordingly. In particular,
1724 the flux density (reported as Integrated Flux in the logger and associated with the "flux" key
1725 in the returned component subdictionary(ies)) for such a case represents the sum of pixel values.
1727 Note also that converting the returned results subdictionary to a component list via cl.fromrecord() currently
1728 only works properly if the flux density units in the results dictionary are conformant with Jy.
1729 If you need to be able to run cl.fromrecord() on the resulting dictionary you can first modify the
1730 flux density units by hand to be (some prefix)Jy and then run cl.fromrecord() on that dictionary,
1731 bearing in mind your unit conversion.
1733 If the input image has units of K, the flux density of components will be reported in units
1734 of [prefix]K*rad*rad, where prefix is an SI prefix used so that the numerical value is between
1735 1 and 1000. To convert to units of K*beam, determine the area of the appropriate beam,
1736 which is given by pi/(4*ln(2))*bmaj*bmin, where bmaj and bmin are the major and minor axes
1737 of the beam, and convert to steradians (=rad*rad). This value is included in the beam portion
1738 of the component subdictionary (key 'beamster'). Then divide the numerical value of the
1739 logged flux density by the beam area in steradians. So, for example
1741 begin{verbatim}
1742 # run on an image with K brightness units
1743 res = imfit(...)
1744 # get the I flux density in K*beam of component 0
1745 comp = res['results']['component0']
1746 flux_density_kbeam = comp['flux']['value'][0]/comp['beam']['beamster']
1747 end{verbatim}
1749 FITTING OVER MULTIPLE CHANNELS
1751 For fitting over multiple channels, the result of the previous successful fit is used as
1752 the estimate for the next channel. The number of gaussians fit cannot be varied on a channel
1753 by channel basis. Thus the variation of source structure should be reasonably smooth in
1754 frequency to produce reliable fit results.
1756 MASK SPECIFICATION
1758 Mask specification can be done using an LEL expression. For example
1760 mask = '"myimage">5' will use only pixels with values greater than 5.
1762 INCLUDING AND EXCLUDING PIXELS
1764 Pixels can be included or excluded from the fit based on their values
1765 using these parameters. Note that specifying both is not permitted and
1766 will cause an error. If specified, both take an array of two numeric
1767 values.
1769 ESTIMATES
1771 Initial estimates of fit parameters may be specified via an estimates
1772 text file. Each line of this file should contain a set of parameters for
1773 a single gaussian. Optionally, some of these parameters can be fixed during
1774 the fit. The format of each line is
1776 peak intensity, peak x-pixel value, peak y-pixel value, major axis, minor axis, position angle, fixed
1778 The fixed parameter is optional. The peak intensity is assumed to be in the
1779 same units as the image pixel values (eg Jy/beam). The peak coordinates are specified
1780 in pixel coordinates. The major and minor axes and the position angle are the convolved
1781 parameters if the image has been convolved with a clean beam and are specified as quantities.
1782 The fixed parameter is optional and is a string. It may contain any combination of the
1783 following characters 'f' (peak intensity), 'x' (peak x position), 'y' (peak y position),
1784 'a' (major axis), 'b' (axial ratio, R = (major axis FWHM)/(minor axis FWHM)),
1785 'p' (position angle). NOTE: One cannot hold the minor axis fixed without holding the major
1786 axis fixed. If the major axis is not fixed, specifying "b" in the fixed string will hold
1787 the axial ratio fixed during the fit.
1789 In addition, lines in the file starting with a # are considered comments.
1791 An example of such a file is:
1793 begin{verbatim}
1794 # peak intensity must be in map units
1795 120, 150, 110, 23.5arcsec, 18.9arcsec, 120deg
1796 90, 60, 200, 46arcsec, 23arcsec, 140deg, fxp
1797 end{verbatim}
1799 This is a file which specifies that two gaussians are to be simultaneously fit,
1800 and for the second gaussian the specified peak intensity, x position, and position angle
1801 are to be held fixed during the fit.
1803 ERROR ESTIMATES
1805 Error estimates are based on the work of Condon 1997, PASP, 109, 166. Key assumptions made are:
1807 - The given model (elliptical Gaussian, or elliptical Gaussian plus constant offset) is an
1808 adequate representation of the data
1809 - An accurate estimate of the pixel noise is provided or can be derived (see above). For the
1810 case of correlated noise (e.g., a CLEAN map), the fit region should contain many "beams" or
1811 an independent value of rms should be provided.
1812 - The signal-to-noise ratio (SNR) or the Gaussian component is large. This is necessary because
1813 a Taylor series is used to linearize the problem. Condon (1997) states that the fractional
1814 bias in the fitted amplitude due to this assumption is of order 1/(S*S), where S is the overall
1815 SNR of the Gaussian with respect to the given data set (defined more precisely below). For a 5
1816 sigma "detection" of the Gaussian, this is a 4% effect.
1817 - All (or practically all) of the flux in the component being fit falls within the selected region.
1818 If a constant offset term is simultaneously fit and not fixed, the region of interest should be
1819 even larger. The derivations of the expressions summarized in this note assume an effectively
1820 infinite region.
1822 Two sets of equations are used to calculate the parameter uncertainties, based on if
1823 the noise is correlated or uncorrelated. The rules governing which set of equations are
1824 used have been described above in the description of the noisefwhm parameter.
1826 In the case of uncorrelated noise, the equations used are
1828 f(A) = f(I) = f(M) = f(m) = k*s(x)/M = k*s(y)/m = (s(p)/sqrt(2))*((M*M - m*m)/(M*m))
1829 = sqrt(2)/S
1831 where s(z) is the uncertainty associated with parameter z, f(z) = s(z)/abs(z) is the
1832 fractional uncertainty associated with parameter z, A is the peak intensity, I is the flux
1833 density, M and m are the FWHM major and minor axes, p is the position angle of the
1834 component, and k = sqrt(8*ln(2)). s(x) and s(y) are the direction
1835 uncertainties of the component measured along the major and minor axes; the resulting
1836 uncertainties measured along the principle axes of the image direction coordinate are
1837 calculated by propagation of errors using the 2D rotation matrix which enacts the rotation through
1838 the position angle plus 90 degrees. S is the overall signal to noise ratio of the component,
1839 which, for the uncorrelated noise case is given by
1841 S = (A/(k*h*r))*sqrt(pi*M*m)
1843 where h is the pixel width of the direction coordinate and r is the rms noise (see the
1844 discussion above for the rules governing how the value of r is determined).
1846 For the correlated noise case, the same equations are used to determine the uncertainties
1847 as in the uncorrelated noise case, except for the uncertainty in I (see below). However,
1848 S is given by
1850 S = (A/(2*r*N)) * sqrt(M*m) * (1 + ((N*N/(M*M)))**(a/2)) * (1 + ((N*N/(m*m)))**(b/2))
1852 where N is the noise-correlation beam FWHM (see discussion of the noisefwhm parameter for
1853 rules governing how this value is determined). "**" indicates exponentiation and a and b
1854 depend on which uncertainty is being calculated. For sigma(A), a = b = 3/2. For M and x,
1855 a = 5/2 and b = 1/2. For m, y, and p, a = 1/2 and b = 5/2. f(I) is calculated in the
1856 correlated noise case according to
1858 f(I) = sqrt( f(A)*f(A) + (N*N/(M*m))*(f(M*f(M) + f(m)*f(m))) )
1860 Note well the following caveats:
1861 - Fixing Gaussian component parameters will tend to cause the parameter uncertainties reported for free
1862 parameters to be overestimated.
1863 - Fitting a zero level offset that is not fixed will tend to cause the reported parameter
1864 uncertainties to be slightly underestimated.
1865 - The parameter uncertainties will be inaccurate at low SNR (a ~10% for SNR = 3).
1866 - If the fitted region is not considerably larger than the largest component that is fit,
1867 parameter uncertainties may be mis-estimated.
1868 - An accurate rms noise measurement, r, for the region in question must be supplied.
1869 Alternatively, a sufficiently large signal-free region must be present in the selected region
1870 (at least about 25 noise beams in area) to auto-derive such an estimate.
1871 - If the image noise is not statistically independent from pixel to pixel, a reasonably accurate noise
1872 correlation scale, N, must be provided. If the noise correlation function is not approximately Gaussian,
1873 the correlation length can be estimated using
1875 N = sqrt(2*ln(2)/pi)* double-integral(dx dy C(x,y))/sqrt(double-integral(dx dy C(x, y) * C(x,y)))
1877 where C(x,y) is the associated noise-smoothing function
1878 - If fitted model components have significan spatial overlap, the parameter uncertainties are likely to
1879 be mis-estimated (i.e., correlations between the parameters of separate components are not accounted
1880 for).
1881 - If the image being analyzed is an interferometric image with poor uv sampling, the parameter
1882 uncertainties may be significantly underestimated.
1884 The deconvolved size and position angle errors are computed by taking the maximum of the absolute values of the
1885 differences of the best fit deconvolved value of the given parameter and the deconvolved size of the eight
1886 possible combinations of (FWHM major axis +/- major axis error), (FWHM minor axis +/- minor axis error),
1887 and (position andle +/- position angle error). If the source cannot be deconvolved from the beam (if the best
1888 fit convolved source size cannot be deconvolved from the beam), upper limits on the deconvolved source size
1889 are sometimes reported. These limits simply come from the maximum major and minor axes of the deconvolved
1890 gaussians taken from trying all eight of the aforementioned combinations. In the case none of these combinations
1891 produces a deconvolved size, no upper limit is reported.
1893 EXAMPLE:
1895 Here is how one might fit two gaussians to multiple channels of a cube using the fit
1896 from the previous channel as the initial estimate for the next. It also illustrates
1897 how one can specify a region in the associated continuum image as the region to use
1898 as the fit for the channel.
1900 begin{verbatim}
1901 imagename = "co_cube.im"
1902 # specify region using region from continuum
1903 region = "continuum.im:source.rgn"
1904 chans = "2~20"
1905 # only use pixels with positive values in the fit
1906 excludepix = [-1e10,0]
1907 # estimates file contains initial parameters for two Gaussians in channel 2
1908 estimates = "initial_estimates.txt"
1909 logfile = "co_fit.log"
1910 # append results to the log file for all the channels
1911 append = "True"
1912 ia.open(imagename)
1913 ia.fitcomponents(region=region, chans=chans, excludepix=excludepix, estimates=estimates, logfile=logfile, append=append)
1914 end{verbatim}
1915 """
1916 return self._swigobj.fitcomponents(box, region, chans, stokes, mask, includepix, excludepix, residual, model, estimates, logfile, append, newestimates, complist, overwrite, dooff, offset, fixoffset, stretch, rms, noisefwhm, summary)
1918 def fromrecord(self, record, outfile=''):
1919 """You can convert an associated image to a record
1920 (torecord) or imagepol tool functions will sometimes give you a record. This function
1921 (fromrecord) allows you to set the contents of an image tool to the content of the record.
1922 This and torecord are used for deserialization and serialization.
1923 """
1924 return self._swigobj.fromrecord(record, outfile)
1926 def getchunk(self, blc=[ int(-1) ], trc=[ int(-1) ], inc=[ int(1) ], axes=[ int(-1) ], list=False, dropdeg=False, getmask=False):
1927 """This function returns the pixels (or optionally the pixel mask) from the
1928 attached image between blc and trc, inclusive. Images with float, complex
1929 float, double, and complex double precision pixel values are supported.
1930 An increment may be specified with inc. Note that if too many pixel values
1931 are retrieved, swapping may occur, result in a decrease in performance,
1932 since the pixel values are stored in memory.
1934 Any illegal blc values are set to zero. Any illegal trc values are set to
1935 the end of the image. If any trc values are less than corresponding blc
1936 values, all the pixel values for that axis are returned. Any illegal inc
1937 values are set to unity.
1939 The axes parameter can be used to reduce the dimensionality of the output
1940 array. It specifies which pixel axes of the image over which to average
1941 the data. For example, consider a 3-D image, with axes=[0,1] and all
1942 other parameters set to their defaults. The result would be a 1-D vector,
1943 a profile along the third axis, with the data averaged over the first two
1944 axes.
1946 A related function is getregion(), which retrieves the pixel values or
1947 pixel mask from a potentially more complex region. Method getchunk() is
1948 retained because it is faster and therefore preferable for repeated
1949 operations in loops if the pixel mask is not required and the region is a
1950 simple box.
1952 If getmask=True, the return value is the pixel mask values, rather than
1953 the pixel values.
1954 """
1955 return self._swigobj.getchunk(blc, trc, inc, axes, list, dropdeg, getmask)
1957 def getregion(self, region={ }, axes=[ int(-1) ], mask='', list=False, dropdeg=False, getmask=False, stretch=False):
1958 """This function recovers the image pixel or pixel mask values in the given region
1959 of interest. Regardless of the shape of the specified, the shape of the pixels and
1960 pixelmask arrays must necessarily be the bounding box of the specified region. If
1961 the region extends beyond the image, it is truncated.
1963 Recall that the recovered pixelmask will reflect both the pixelmask stored in the
1964 image, and the region (their masks are 'anded' together).
1966 The argument axes can be used to reduce the dimensionality of the output array. It
1967 specifies which pixel axes of the image to average the data over. For example,
1968 consider a 3-D image. With axes=[0,1] and all other arguments left at their
1969 defaults, the result would be a 1-D vector, a profile along the third axis, with
1970 the data averaged over the first two axes.
1972 This method differs in a couple of ways from the getchunk() method. For example,
1973 the specified region can be much more complex (eg, a union of polygons) than the
1974 limited, simple regions that can be specified in getchunk(), which must be
1975 rectangular. On the other hand, getregion() is less effective than getchunk()
1976 for the same region specification. So if one is interested in say, iterating
1977 through an image, getting a regular hyper-cube of pixels and doing something
1978 with them, then getchunk() will be faster. This would be especially noticeable if
1979 you iterated line by line through a large image (and of course, in both cases,
1980 retrieving very large regions will become very resource intensive, as these
1981 returned arrays are completely stored in memory).
1982 """
1983 return self._swigobj.getregion(region, axes, mask, list, dropdeg, getmask, stretch)
1985 def getprofile(self, axis=int(-1), function='mean', region={ }, mask='', unit='', stretch=False, spectype='default', restfreq='', frame='', logfile=''):
1986 """This application returns information on a one-dimensional profile taken along a specified image axis.
1987 The region of interest is collapsed (a'la ia.collapse() along all axes orthogonal to the one specified, and)
1988 the specified aggregate function is applied to these pixels to generate the returned values.
1990 The aggregate function must be one of the functions supported by ia.collapse; ie, 'flux', 'madm', 'max', 'mean',
1991 'median', 'min', 'rms', 'stdev', 'sum', 'variance', and 'xmadm'. See the help for ia.collapse() for details regarding
1992 these functions. Minimum match and case insenstivity is supported. In addition, single binary (addition,
1993 subtraction, multiplication, and division) operations of these functions are supported, eg function='max*min'
1994 will return data that is the product of the maximum and the mininum for each plane along the specified
1995 axis.
1997 One may specify the unit of the returned coordinate values. Unless axis is the spectral axis, unit must be
1998 conformant with the corresponding axis unit in the image coordinate system or it must be 'pixel' which signifies,
1999 pixel, rather than world, coordinate values should be calculated. If axis is the spectral axis, unit may be a
2000 velocity unit (assuming the coordinate system has a rest frequency or restfreq is specified) or a length unit.
2001 In these cases, the returned coordinate values will be converted to velocity or wavelength, respectively.
2003 The parameter spectype may be used to specify the velocity or wavelength type for the returned coordinate values
2004 if profile is taken along spectral axis. Supported (minimum match, case insensitive) values) are "relativistic
2005 velocity", "beta", "radio velocity", "optical velocity", "wavelength", "air wavelength", "default". The "default"
2006 value is equivalent to "relativistic" if unit is a velocity unit or "wavelength" if unit is a length unit.
2008 The restfreq parameter allows one to set the rest frequency for the coordinates to be returned if axis is the
2009 spectral axis and unit is a velocity unit. If blank, the rest frequency associated with the image coordinate
2010 system is used.
2012 The frame allows one to specify which kinematic reference frame that the returned coordinate values should be
2013 calculated in. It is only used if axis is the spectral axis and unit is unspecified or is specified and a
2014 frequency unit. If blank, the reference frame associated with the image coordinate system is used.
2016 The returned dictionary
2017 contains the keys:
2019 values: one-dimensional array along the specified axis containing values resulting from applying the specified
2020 aggregate function to corresponding pixels at the same location along that axis.
2021 mask: one-dimensional array of booleans of the resulting mask after applying the aggregate function, formed in the
2022 same way as that formed by ia.collapse.
2023 coords One-dimensional array of corresponding coordinate values along the specified axis in the specified unit
2024 (or the unit associated with the axis in the image coordinate system if unspecified).
2025 xUnit The unit used for calculating the values the coords array.
2026 """
2027 return self._swigobj.getprofile(axis, function, region, mask, unit, stretch, spectype, restfreq, frame, logfile)
2029 def getslice(self, x, y, axes=[ int(0),int(1) ], coord=[ int(-1) ], npts=int(0), method='linear'):
2030 """This function returns a 1-D slice (the pixels and opionally the pixel mask) from
2031 the attached image. The slice is constrained to lie in a plane of two cardinal
2032 axes (e.g. XY or YZ). Interpolation is permitted between pixels, and a set of
2033 interpolation schemes is available.
2035 The slice is specified as a polyline giving the x and y coordinates and the axes
2036 of the plane holding that slice. The absolute pixel coordinates of the other
2037 axes may be specified using the coord parameter. If not specified, these values
2038 default to pixel 0 on the relevant axes.
2040 The npts parameter allows the number of values to be returned.
2042 The method parameter allows specification of the interpolation method to be
2043 used. Allowed values are 'nearest', 'linear', and 'cubic'. In the case of an
2044 image with complex valued pixels, the interpolation is done independently on the
2045 real and imaginary values. For example, the linearly interpolated value midway
2046 between pixels with values of 1 + 5j and 2 + 7j would be 1.5 + 6j.
2048 The return value is a dictionary with keys 'pixels' (interpolated pixel values),
2049 'mask' (interpolated mask), 'xpos' (x-location in absolute pixel coordinates),
2050 'ypos' (y-location in absolute pixel coordinates), 'distance' (distance along
2051 slice in pixels), and 'axes' (the x and y axes of slice).
2052 """
2053 return self._swigobj.getslice(x, y, axes, coord, npts, method)
2055 def hanning(self, outfile='', region={ }, mask='', axis=int(-10), drop=True, overwrite=False, stretch=False, dmethod='copy'):
2056 """This application performs Hanning convolution of one axis of an image defined by
2058 z[i] = 0.25*y[i-1] + 0.5*y[i] + 0.25*y[i+1] (equation 1)
2060 where z[i] is the value at pixel i in the hanning smoothed image, and
2061 y[i-1], y[i], and y[i+1] are the values of the input image at pixels i-1,
2062 i, and i+1 respectively. It supports both float and complex valued images.
2063 The length of the axis along which the convolution is to occur must be at least
2064 three pixels in the selected region. Masked pixel values are set to zero prior to
2065 convolution. All nondefault pixel masks are ignored during the calculation.
2067 The convolution is done in the image domain (i.e., not with an FFT).
2069 If drop=False, the length of the output axis will be the same as that of the input
2070 axis. The output pixel values along the convolution axis will be related to those
2071 of the input values according to equation 1, except the first and last pixels. In that
2072 case,
2074 z[0] = 0.5*(y[0] + y[1])
2076 and,
2078 z[N-1] = 0.5*(y[N-2] + y[N-1])
2080 where N is the number of pixels along the convolution aixs.
2081 The pixel mask, ORed with the OTF mask if specified, is copied from the selected
2082 region of the input image to the output image. Thus for example, if the selected
2083 region in the input image has six planes along the convolution axis, and if the pixel
2084 values, which are all unmasked, on a slice along this axis are [1, 2, 5, 10, 17, 26],
2085 the corresponding output pixel values will be [1.5, 2.5, 5.5, 10.5, 17.5, 21.5].
2087 If drop=True and dmethod="copy", the output image is the image calculated if
2088 drop=True, except that only the odd-numbered planes are kept. Furthermore, if the
2089 number of planes along the convolution axis in the selected region of the input image
2090 is even, the last odd number plane is also discarded. Thus, if the selected region
2091 has N pixels along the convolution axis in the input image, along the convolution
2092 axis the output image will have (N-1)/2 planes if N is odd, or (N-2)/2 planes if N
2093 is even. In this case, the pixel and mask values are copied directly, without further
2094 processing. Thus for example, if the selected region in the input image has six planes
2095 along the convolution axis, and if the pixel values, which are all unmasked, on a slice
2096 along this axis are [1, 2, 5, 10, 17, 26], the corresponding output pixel values will be
2097 [2.5, 10.5].
2099 If drop=True and dmethod="mean", first the image described in the drop=False case
2100 is calculated. The first plane and last plane(s) of that image are then discarded as
2101 described in the drop=True, dmethod="copy" case. Then, the ith plane of the output
2102 image is calculated by averaging the (2*i)th and (2*i + 1)th planes of the intermediate
2103 image.Thus for example, if the selected region in the input image has six planes
2104 along the convolution axis, and if the pixel values, which are all unmasked, on a slice
2105 along this axis are [1, 2, 5, 10, 17, 26], the corresponding output pixel values will be
2106 [4.0, 14.0]. Masked values are taken into consideration when forming this average, so if
2107 one of the values is masked, it is not used in the average. If at least one of the values
2108 in the input pair is not masked, the corresponding output pixel will not be masked.
2110 The hanning smoothed image is written to disk with name {stfaf outfile}, if
2111 specified. If not, no image is written but the image is still accessible via
2112 the returned image analysis tool (see below).
2114 This method always returns an image analysis tool which is attached to the
2115 hanning smoothed image. This tool should always be captured and closed after
2116 any desired manipulations have been done. Closing the tool frees up system
2117 resources (eg memory), eg,
2119 hanning_image = ia.hanning(...)
2120 begin{verbatim}
2121 # do things (or not) with hanning_image
2122 ...
2123 # close the returned tool promptly upon finishing with it.
2124 end{verbatim}
2125 hanning_image.done()
2127 See also the other convolution functions
2128 convolve2d,
2129 sepconvolve and
2130 convolve.
2131 """
2132 return _wrap_image(swig_object=self._swigobj.hanning(outfile, region, mask, axis, drop, overwrite, stretch, dmethod))
2134 def haslock(self):
2135 """This function can be used to find out whether the image has a read or a
2136 write lock set. It is not of general user interest. It returns
2137 a vector of Booleans of length 2. Position 1 says whether
2138 a read lock is set, position 2 says whether a write lock is set.
2140 In general locking is handled automatically, with a built in lock
2141 release cycle. However, this function can be useful in scripts when a
2142 file is being shared between more than one process. See also functions
2143 unlock and
2144 lock.
2145 """
2146 return self._swigobj.haslock()
2148 def histograms(self, axes=[ int(-1) ], region={ }, mask='', nbins=int(25), includepix=[ float(-1) ], cumu=False, log=False, stretch=False):
2149 """This method computes histograms of the pixel values in the image.
2150 The values are returned in a dictionary.
2152 The chunk of the image over which you compute the histograms is
2153 specified by a vector of axis numbers (argument {stfaf axes}). For
2154 example, consider a 3-dimensional image for which you specify {stfaf
2155 axes=[0,2]}. The histograms would be computed for each XZ (axes 0 and
2156 2) plane in the image. You could then examine those histograms as a
2157 function of the Y (axis 1) axis. Or perhaps you set {stfaf axes=[2]},
2158 whereupon you could examine the histogram for each Z (axis 2) profile as
2159 a function of X and Y location in the image.
2161 You have control over the number of bins for each histogram ({stfaf
2162 nbins}). The bin width is worked out automatically for each histogram
2163 and may vary from histogram to histogram (the range of pixel values is
2164 worked out for each chunk being histogrammed).
2166 You have control over which pixels are included in the histograms via
2167 the {stfaf includepix} argument. This vector specifies a range of
2168 pixel values to be included in the histograms. If you only give one
2169 value for this, say {stfaf includepix=[b]}, then this is interpreted as
2170 {stfaf includepix=[-abs(b),abs(b)]}. If you specify an inclusion
2171 range, then the range of pixel intensities over which the histograms are
2172 binned is given by this range too. This is a way to make the bin width
2173 the same for each histogram.
2175 You can control if the histogram is cumulative or non-cumulative via the
2176 cumu parameter.
2178 You have countrol over how the bin counts are returned. If log = false,
2179 the actual counts are returned. If true, the values returned are the log10
2180 values of the actual counts.
2182 The results are returned as a dictionary. The counts (field "counts") and
2183 the abscissa values (field "values") for all bins in each histogram are returned.
2184 The shape of the first dimension of those arrays contained in those fields is {stfaf nbins}.
2185 The number and shape of the remaining dimensions are those of the display axes(the
2186 axes in the image for which you did not compute the histograms). For example, if one
2187 has a three dimensional image and sets {stfaf axes=[2]}, the display axes are 0 and 1,
2188 so the shape of each counts and values array is then [nbins,nx,ny], where nx and ny
2189 are the length of the zeroth and first axes, respectively.
2191 In addition, the mean (field "mean") and standard deviation (field "sigma") computed
2192 using the data in each histogram is returned. The shape of these arrays is equal to
2193 the shape of the display axes. So,
2194 """
2195 return self._swigobj.histograms(axes, region, mask, nbins, includepix, cumu, log, stretch)
2197 def history(self, list=True):
2198 """This method interogates the history of an image.
2200 The history is returned as an array of strings, where each element represents
2201 an individual history entry. If True, the list parameter will also cause the
2202 history to be emitted by the logger.
2204 Note that entries can be permanently added to the image history by using the
2205 ia.sethistory() method.
2206 """
2207 return self._swigobj.history(list)
2209 def insert(self, infile='', region={ }, locate=[ float(-1) ], verbose=False):
2210 """This function inserts the specified image (or part of it) into the image
2211 referenced by this tool.
2212 The specified image may be given via argument {stfaf infile}
2213 as a disk file name (it may be in native casa, fits, or Miriad
2214 format; Look htmlref{here}{IMAGES:FOREIGNIMAGES} for more
2215 information on foreign images).
2217 If the {stfaf locate} vector is not given, then the images are
2218 aligned (to an integer pixel shift) by their reference pixels.
2220 If {stfaf locate} vector is given, then those values that are given,
2221 give the absolute pixel in the output (this) image of the bottom left
2222 corner of the input (sub)image. For those values that are not given,
2223 the input image is symmetrically placed in the output image.
2225 The image referenced by this tool is modified in place; no new image
2226 is created. The method returns True if successful.
2227 """
2228 return self._swigobj.insert(infile, region, locate, verbose)
2230 def isopen(self):
2231 """This method returns True if the image tool has an attached image.
2232 """
2233 return self._swigobj.isopen()
2235 def ispersistent(self):
2236 """This function can be used to find out whether the image is persistent on
2237 disk or not. There is a subtle difference from the image being
2238 virtual. For example, a virtual image which references another
2239 which is on disk is termed persistent.
2240 """
2241 return self._swigobj.ispersistent()
2243 def lock(self, writelock=False, nattempts=int(0)):
2244 """This function can be used to acquire a Read or a Read/Write lock
2245 on the imagefile. It is not of general user interest.
2247 In general locking is handled automatically, with a built in lock
2248 release cycle. However, this function can be useful in scripts when a
2249 file is being shared between more than one process. See also functions
2250 unlock and haslock.
2251 """
2252 return self._swigobj.lock(writelock, nattempts)
2254 def makecomplex(self, outfile, imag, region={ }, overwrite=False):
2255 """This function combines the current image with another image to make
2256 a complex image. The current image (i.e. that associated with this
2257 Image must have real valued pixels). The image used for generating the
2258 imaginary part of the pixel values is specified using the imag parameter, and
2259 it must persistent. The image attached to this tool and the image specified
2260 using the imag parameter must have the same precision, or else an exception
2261 will be thrown. If both are float precision, the resulting image will have
2262 float precision pixel values. If both are double precision, the resulting image
2263 will have double precision pixel values. The coordinate systems of the two
2264 input images must be conformant. The metadata written to the resulting image is
2265 copied from the image attached to this tool.
2266 """
2267 return self._swigobj.makecomplex(outfile, imag, region, overwrite)
2269 def maskhandler(self, op='default', name=[ ]):
2270 """This method is used to manage or handle pixel masks . A CASA image may contain
2271 zero, one, or more pixel masks. Any of these masks can be designated the
2272 default pixel mask assoicated with the given image. The default mask is acted
2273 upon and/or used by CASA applications. For example, if ia.statistics() will
2274 exclude pixels which are masked as bad (False) from the calculations.
2276 This method does not modify the individual boolean values of any masks.
2278 The op parameter is used to specify the behaviour. In all cases, the specified
2279 operation can be specified by a three character string. Supported values of the
2280 op parameter are:
2282 'default': Retrieve the name of the default pixel mask associated with the
2283 image. A one element array containing the empty string is returned if the image
2284 has no default mask.
2286 'get': Retrieve the name(s) of all the existing pixel masks. Note that the
2287 ia.summary() method may also be used to view the pixel masks associated with an
2288 image.
2290 'set': Set the default pixel mask to the value specified by the name parameter.
2291 If the value of the name parameter is the empty string, then the default mask is
2292 unset (ie, all the pixels will be treated as being unmasked).
2294 'delete': Delete the pixel mask(s) specified by the name parameter. To delete
2295 more than one mask, the name parameter can be an array of strings. Any supplied
2296 pixel mask name that does not exist is silently ignored.
2298 'rename': Rename the mask specified as the first element of the name array
2299 parameter (name[0]) to the value specified in the second element of the name
2300 array parameter (name[1]]. In this case, the name array parameter must have
2301 exactly two elements.
2303 'copy': Copy the mask specified in the first element of the name array
2304 parameter (name[0]) to the a mask name specified in the second element of the
2305 name array parameter (name[1]]. In this case, the name array parameter must have
2306 exactly two elements. A mask from another image can be copied by using the
2307 imagename:maskname syntax for the first element in the name array, eg,
2308 'myimage:mask2'.
2309 """
2310 return self._swigobj.maskhandler(op, name)
2312 def miscinfo(self):
2313 """A casa imagefile can accumulate miscellaneous information
2314 during its lifetime. This information is stored in a record called the {stff
2315 miscinfo} record. For example, the fits filler puts header keywords
2316 it doesn't otherwise use into the {stff miscinfo} record. This {stff
2317 miscinfo} record is not guaranteed to have any entries, so it's up to
2318 you to check for any fields that you require.
2320 You can also put things into this record (see
2321 setmiscinfo) yourself, to keep
2322 information that the system might not otherwise store for you.
2324 When the image is written out to fits, the items in the
2325 {stff miscinfo} record are written to the fits file
2326 as keywords with the corresponding record field name.
2327 """
2328 return self._swigobj.miscinfo()
2330 def modify(self, model, region={ }, mask='', subtract=True, list=True, stretch=False):
2331 """This function applies a model of the sky to the image. You can add or
2332 subtract the model which is contained in a
2333 Componentlist tool.
2335 The pixel values are only changed where the total mask
2336 (combination of the default pixel mask [if any] and the OTF mask)
2337 is good (True). If the computation fails for a particular
2338 pixel (e.g. coordinate undefined) that pixel will be
2339 masked bad.
2341 DISK MODELS
2343 Pixels with centers inside the disk will have the same values, even if a pixel straddles the
2344 edge of the disk. Pixels with straddle the edge of the disk which have centers outside the
2345 disk are given values of zero. Thus, one should not expect the flux density of the disk to
2346 be exactly the provided value to the component list; for a given size disk, the computed flux
2347 density will be closer to the expected value for images with smaller pixels.
2348 """
2349 return self._swigobj.modify(model, region, mask, subtract, list, stretch)
2351 def maxfit(self, region={ }, point=True, width=int(5), negfind=False, list=True):
2352 """This function finds the pixel with the maximum value in the region, and
2353 then uses function findsources
2354 to generate a Componentlist with one component. The component
2355 will be of type Point ({stfaf point=T}) or Gaussian ({stfaf point=F}).
2357 If {stfaf negfind=F} the maximum pixel value is found in the region and fit.
2358 If {stfaf negfind=T} the absolute maximum pixel value is found in the region
2359 and fit.
2361 See function findsources for
2362 a description of arguments {stfaf point} and {stfaf width}.
2364 See also the function fitcomponents.
2365 """
2366 return self._swigobj.maxfit(region, point, width, negfind, list)
2368 def moments(self, moments=[ int(0) ], axis=int(-10), region={ }, mask='', method=[ ], smoothaxes=[ int(-1) ], smoothtypes=[ ], smoothwidths=[ float(0.0) ], includepix=[ float(-1) ], excludepix=[ float(-1) ], peaksnr=float(3.0), stddev=float(0.0), doppler='RADIO', outfile='', smoothout='', overwrite=False, drop=True, stretch=False):
2369 """noindent{bf Summary}
2371 The primary goal of this function is to enable you to analyze a
2372 multi-dimensional image by generating moments of a specified axis.
2373 This is a time-honoured spectral-line analysis technique used for
2374 extracting information about spectral lines.
2376 You can generate one or more output moment images. The return value
2377 of this function is an on-the-fly Image tool holding the {bf first}
2378 of the output moment images.
2380 The word 'moment' is used loosely here. It refers to collapsing an axis
2381 (the moment axis) to one pixel and setting the value of that pixel (for
2382 all of the other non-collapsed axes) to something computed from the data
2383 values along the moment axis. For example, take an RA-DEC-Velocity
2384 cube, collapse the velocity axis by computing the mean intensity at each
2385 RA-DEC pixel. This function offers many different moments and a variety
2386 of automatic methods to compute them.
2388 We try to make a distinction between a 'moment' and a 'method'. This
2389 boundary is a little blurred, but it claims to refer to the distinction
2390 between what you are computing, and how the pixels that were included in
2391 that computation were selected. For example, a 'moment' would be the
2392 average value of some pixel values in a spectrum. A 'method' for
2393 selecting those pixels would be a simple pixel value range specifying
2394 which pixels should be included.
2396 There are many available moments, and you specify each one with an
2397 integer code as it would get rather cumbersome to refer to them via
2398 strings. In the list below, the value of the $i$th pixel of the
2399 spectrum is $I_i$, the coordinate of this pixel is $v_i$ (of course it
2400 may not be velocity), and there are $n$ pixels in the spectrum. The
2401 available moments are:
2403 begin{itemize}
2404 item{$-1$} -- the mean value of the spectrum
2405 begin{displaymath}
2406 { {1over n} {sum {I_i}}}
2407 end{displaymath}
2408 medskip
2410 item{0} -- the integrated value of the spectrum
2411 begin{displaymath}
2412 M_0 = Delta v sum I_i
2413 end{displaymath}
2415 where $Delta v$ is the width (in world coordinate units) of a pixel
2416 along the moment axis
2417 medskip
2419 item{1} -- the intensity weighted coordinate (this is
2420 traditionally used to get 'velocity fields')
2422 begin{displaymath}
2423 M_1 = { {sum {I_i v_i}} over {M_0}}
2424 end{displaymath}
2425 medskip
2427 item{2} -- the intensity weighted dispersion of the coordinate
2428 (this is traditionally used to get 'velocity dispersion fields')
2430 begin{displaymath}
2431 sqrt{ { {sum {I_i left(v_i - M_1right)^2}} over {M_0}}}
2432 end{displaymath}
2433 medskip
2435 item{3} -- the median of $I$
2436 medskip
2438 item{4} -- the median coordinate. Here we treat the spectrum as a
2439 probability distribution, generate the cumulative distribution, and then
2440 find the coordinate corresponding to the 50% value. This moment is not
2441 very robust, but it is useful for quickly generating a velocity field in
2442 a way that is not sensitive to noise. However, it will only give
2443 sensible results under certain conditions. The generation of the
2444 cumulative distribution and the finding of the 50% level really only
2445 makes sense if the cumulative distribution is monotonic. This
2446 essentially means only selecting pixels which are positive or negative.
2447 For this reason, this moment type is only supported with the basic
2448 method (see below -- i.e. no smoothing, no windowing, no fitting) with
2449 a pixel selection range that is either all positive, or all negative
2450 medskip
2452 item{5} -- the standard deviation about the mean of the spectrum
2453 begin{displaymath}
2454 sqrt{ {1over {left(n-1right)}} sum{left(I_i - bar{I}right)^2 }}
2455 end{displaymath}
2456 medskip
2458 item{6} -- the root mean square of the spectrum
2459 begin{displaymath}
2460 sqrt{ {1 over n} sum{I_i^2}}
2461 end{displaymath}
2462 medskip
2464 item{7} -- the absolute mean deviation of the spectrum
2465 begin{displaymath}
2466 {1 over n} sum {|(I_i - bar{I})|}
2467 end{displaymath}
2468 medskip
2470 item{8} -- the maximum value of the spectrum
2471 medskip
2472 item{9} -- the coordinate of the maximum value of the spectrum
2473 medskip
2474 item{10} -- the minimum value of the spectrum
2475 medskip
2476 item{11} -- the coordinate of the minimum value of the spectrum
2477 medskip
2478 end{itemize}
2480 bigskip
2481 noindent {Smoothing}
2483 The purpose of the smoothing functionality is purely to provide
2484 a mask. Thus, you can smooth the input image, apply a pixel
2485 include or exclude range, and generate a smoothed mask which is then
2486 applied before the moments are generated. The smoothed data
2487 are not used to compute the actual moments; that is always done
2488 from the original data.
2490 bigskip
2491 noindent{bf Basic Method}
2493 The basic method is to just compute moments directly from the pixel
2494 values. This can be modified by applying pixel value inclusion or
2495 exclusion ranges (arguments {stfaf includepix} and {stfaf excludepix}).
2497 You can then also convolve the image (arguments {stfaf smoothaxes}, {stfaf
2498 smoothtypes}, and {stfaf smoothwidths}) and find a mask based on the inclusion
2499 or exclusion ranges applied to the convolved image. This mask is then
2500 applied to the unsmoothed data for moment computation.
2502 bigskip
2503 noindent{bf Window Method}
2505 The window method (invoked with argument {stfaf method='window'}) does
2506 no pixel-value-based selection. Instead a window is found (hopefully
2507 surrounding the spectral line feature) and only the pixels in that
2508 window are used for computation. This window can be found from the
2509 convolved or unconvolved image (arguments {stfaf smoothaxes}, {stfaf
2510 smoothtypes}, and {stfaf smoothwidths}).
2512 The moments are always computed from the unconvolved data. The window
2513 can be found (for each spectrum) automatically. The
2514 automatic methods are via Bosma's converging mean algorithm ({stfaf
2515 method='window'}) or by fitting Gaussians and taking $pm 3sigma$ as
2516 the window ({stfaf method='window,fit'}).
2517 In Bosma's algorithm, an initial guess for a range of pixels surrounding
2518 a spectral feature is refined by widening until the mean of the pixels
2519 outside of the range converges (to the noise).
2521 bigskip
2522 noindent{bf Fit Method}
2524 The fit method ({stfaf method='fit'}) fits Gaussians to spectral
2525 features automatically. The moments are then computed from the
2526 Gaussian fits (not the data themselves).
2528 bigskip
2529 noindent{bf Other Arguments}
2531 begin{itemize}
2533 item {stfaf outfile} - If you are creating just one moment image,
2534 and you specify {stfaf outfile}, then the image is created
2535 on disk with this name. If you leave {stfaf outfile} empty
2536 then a temporary image is created. In both cases, you can
2537 access this image with the returned Image tool. If you are
2538 making more than one moment image, then theses images are always
2539 created on disk. If you specify {stfaf outfile} then this is
2540 the root for the output file names. If you don't specify it,
2541 then the input image name is used as the root.
2543 item {stfaf smoothing} - If you smooth the image to generate a
2544 mask, you specify the kernel widths via the {stfaf smoothwidths}
2545 argument in the same way as in the
2546 sepconvolve function. See it for
2547 details.
2549 item {stfaf stddev} - Some of the automatic methods also require an
2550 estimate of the noise level in the image. This is used to assess
2551 whether a spectrum is purely noise or not, and whether there is any
2552 signal worth digging out. If you don't give it via the {stfaf stddev}
2553 argument, it will be worked out automatically from a Gaussian fit to the
2554 bins above 25% from a histogram of the entire image.
2556 item {stfaf includepix, excludepix} - The vectors given by arguments
2557 {stfaf includepix} and {stfaf excludepix} specify a range of pixel
2558 values for which pixels are either included or excluded. They are
2559 mutually exclusive; you can specify one or the other, but not both. If
2560 you only give one value for either of these, say {stfaf includepix=b},
2561 then this is interpreted as {stfaf includepix=[-abs(b),abs(b)]}.
2563 The convolving point-spread function is normalized to have a volume of
2564 unity. This means that point sources are depressed in value, but
2565 extended sources that are large with respect to the PSF remain
2566 essentially on the same intensity scale; these are the structures you
2567 are trying to find with the convolution so this is what you want.
2568 If you convolve the image, then arguments like {stfaf includepix} select
2569 based upon the convolved image pixel values. If you are having trouble
2570 getting these right, you can output the convolved image ({stfaf smoothout})
2571 and assess the validity of your pixel ranges. Note also that if you are
2572 Hanning convolving (usually used on a velocity axis), then the width for
2573 this kernel must be 3 pixels (triangular smoothing kernels of other
2574 widths have no valid theoretical basis).
2576 item {stfaf doppler} - If you compute the moments along a spectral
2577 axis, it is conventional to compute the world coordinate (needed for
2578 moments 0, 1 and 2) along that axis in "km/s". The argument {stfaf
2579 doppler} lets you specify what doppler convention the velocity will be
2580 calculated in. You can choose from {stfaf doppler=radio, optical,
2581 true}. See function summary for the
2582 definitions of these codes. For other moment-axis types, the world coordinate
2583 is computed in the native units.
2585 item {stfaf mask} - The total input mask is the combination of the
2586 default pixelmask (if any) and the OTF mask. Once this mask
2587 has been established, then the moment method may make additional
2588 pixel selections.
2590 item {stfaf drop} - If this is true (the default) then the moment axis
2591 is dropped from the output image. Otherwise, the output images have a
2592 moment axis of unit length and coordinate information that is the same
2593 as for the input image. This coordinate information may be totally
2594 meaningless for the moment images.
2596 end{itemize}
2598 Finally, if you ask for a moment which requires the coordinate to be
2599 computed for each profile pixel (these are the intensity weighted mean
2600 coordinate [moment 1] and the intensity weighted dispersion of the
2601 coordinate [moment 2]), and the profile axis is not separable then there
2602 will be a performance loss. Examples of non-separable axes are RA and
2603 Dec. If the axis is separable (e.g. a spectral axis) there is no
2604 penalty. In the latter case, the vector of coordinates for one profile
2605 is the same as the vector for another profile, and it can be precomputed
2606 (once).
2608 Note that this function has no ``virtual'' output file capability. All
2609 output files are written to disk. The output mask for these images is
2610 good (T) unless the moment method fails to generate a value (e.g. the
2611 total input pixel mask was all bad for the profile) in which case it will be bad (F).
2613 If an image has multiple (per-channel beams) and the moment axis is equal to the
2614 spectral axis, each channel will be convolved with a beam that is equal to the beam
2615 having the largest area in the beamset prior to moment determination.
2616 """
2617 return _wrap_image(swig_object=self._swigobj.moments(moments, axis, region, mask, method, smoothaxes, smoothtypes, smoothwidths, includepix, excludepix, peaksnr, stddev, doppler, outfile, smoothout, overwrite, drop, stretch))
2619 def name(self, strippath=False):
2620 """This function returns the name of the imagefile By default, this
2621 function returns the full absolute path of the imagefile. You can
2622 strip this path off if you wish with the {stfaf strippath} argument and
2623 just recover the imagefile name itself.
2624 """
2625 return self._swigobj.name(strippath)
2627 def open(self, infile, cache=True):
2628 """This method detaches from the current image (if an image is attached to the tool), and
2629 reattaches it (opens) to the new image.
2631 The input image file may be in native CASA, FITS, or Miriad format. In the case
2632 of CASA images, images with float, complex float, double, and complex double
2633 valued pixels are supported. Note that only FITS images with float valued pixels
2634 are supported.
2636 The cache parameter applies only to component list images and indicates if pixel
2637 values should be cached after they are computed for faster retrieval. It is not
2638 used for other image types.
2639 """
2640 return self._swigobj.open(infile, cache)
2642 def pad(self, outfile='', npixels=int(1), value=float(0), padmask=False, overwrite=False, region={ }, box='', chans='', stokes='', mask='', stretch=False, wantreturn=True):
2643 """This method pads the directional plane of an image with a specified number of pixels on each side. The
2644 numerical and mask values of the padding pixels may also be specified. If a region is selected, a subimage
2645 of that region is created and then padded with the specified pixel parameters. Thus, padding an image of
2646 shape (ra, dec, freq) = (512, 512, 10) specifying npixels = 3 results in an image of size (518, 518, 10), with
2647 the blc of the directional plane of the original pixel set corresponding to the directional pixel of (3, 3)
2648 in the output.
2649 If wantreturn is True, an image analysis tool attached to the output image is returned. If False, none is
2650 returned.
2652 """
2653 return _wrap_image(swig_object=self._swigobj.pad(outfile, npixels, value, padmask, overwrite, region, box, chans, stokes, mask, stretch, wantreturn))
2655 def crop(self, outfile='', axes=[ ], overwrite=False, region={ }, box='', chans='', stokes='', mask='', stretch=False, wantreturn=True):
2656 """This method crops masked slices from the perimeter of an image. The axes parameter specifies which axes to
2657 consider. Axes not specified will not be cropped. An empty array implies that all axes should be considered.
2658 If wantreturn is True, an image analysis tool attached to the output image is returned. If False, none is
2659 returned.
2661 """
2662 return _wrap_image(swig_object=self._swigobj.crop(outfile, axes, overwrite, region, box, chans, stokes, mask, stretch, wantreturn))
2664 def pixelvalue(self, pixel=[ int(-1) ]):
2665 """This function gets the value of the image and the mask at the specified
2666 pixel coordinate. The values are returned in a record with fields
2667 'value', 'mask' and 'pixel'. The value is returned as a quantity, the mask
2668 as a Bool (T is good). The 'pixel' field holds the actual
2669 pixel coordinate used.
2671 If the specified pixel coordinate is off the image, "{}" is returned.
2673 Excessive elements in {stfaf pixel} are silently discarded.
2674 Missing elements are given the (nearest integer) value of the reference pixel.
2675 This is reflected in the output record 'pixel' field.
2676 """
2677 return self._swigobj.pixelvalue(pixel)
2679 def putchunk(self, pixels=[ ], blc=[ int(-1) ], inc=[ int(1) ], list=False, locking=True, replicate=False):
2680 """This method puts an array into the image (which must be writable, eg it
2681 will fail on FITS images). If there is a default pixel mask, it is
2682 ignored. It is the complement of the getchunk() method. The blc, trc,
2683 and increment (inc) may be specified. If they are unspecified, they
2684 default to the beginning of the image and an increment of one.
2686 Any illegal blc values are set to zero. Any illegal inc values are set
2687 to unity.
2689 An error will result if an attempt is made to put an array the extends
2690 beyond the image edge (i.e., it is not truncated or decimated).
2692 If there are fewer axes in the array than in the image, the array is
2693 assumed to have trailing axes of length unity. Thus, if you have a 2D
2694 array and want to put it in as the YZ plane rather than the XY plane,
2695 you must ensure that the shape of the array is [1,nx,ny].
2697 However, the replicate parameter can be used to replicate the array
2698 throughout the image (from the blc to the trc). For example, if a 2D
2699 array is provided for a 3D image, it can be replicated along the third
2700 axis by setting replicate=True. The replication is done from the
2701 specified blc to the end of the image. Method putregion() can be used
2702 to terminate the replication at a trc value.
2704 The locking parameter controls two aspects. If True, then after the
2705 method is called, the image is unlocked (so some other process can
2706 acquire a lock) and it is indicated that the image has
2707 changed. The reason for having this argument is that the unlocking and
2708 updating processes are quite expensive. If putchunk() is called
2709 repeatedly in eg, a loop, it is advisable to set this parameter to True.
2711 A related function is putregion(), which supports putting the pixel and
2712 mask values into a more complex region. Method putchunk() is faster and
2713 therefore preferable for repeated operation in loops if the pixel mask
2714 is not required.
2716 See also the methods set() and calc() which can also be used to change
2717 pixel values.
2718 """
2719 return self._swigobj.putchunk(pixels, blc, inc, list, locking, replicate)
2721 def putregion(self, pixels=[ ], pixelmask=[ ], region={ }, list=False, usemask=True, locking=True, replicate=False):
2722 """This function replaces data and/or pixel mask values in the image in the
2723 specified region. The pixels and/or pixelmask arrays must be the shape of
2724 the bounding box, and the whole bounding box is replaced in the image. The
2725 region is only used to specify the bounding box. If the region extends
2726 beyond the image, it is truncated. If the pixels or pixelmask array shapes
2727 do not match the bounding box, an error will result. Values in the pixels
2728 array must share the same domain as the pixel values in the image. If the
2729 pixels array contains real values and the image pixels contain complex
2730 values (or vice versa), an exception will be thrown.
2732 When you put a pixel mask, it either replaces the current default pixel mask,
2733 or is created.
2735 The usemask parameter is only relevant when you are putting pixel values and
2736 there is a pixel mask (meaning also the one you might have just put in place).
2737 If usemask=True, then only pixels for which the mask is good (True) are
2738 altered. If usemask=False, then all the pixels in the region are altered.
2740 The replicate parameter can be used to replicate the array throughout the
2741 image (from the blc to the trc). For example, if a two dimensional array is
2742 provided for a three dimensional image, it can be replicated along the third
2743 axis by setting replicate=True. The replication is done in the specified
2744 region.
2746 The locking parameter controls two things. If True, then after the method
2747 is called, the image is unlocked (so some other process can acquire a lock)
2748 and it is indicated that the image has changed. The reason for this
2749 parameter is that the unlocking and updating processes are quite expensive.
2750 If putregion() is being called multiple times, in a for loop, for example,
2751 it is recommended to set locking=True (and to consider using putchunk()
2752 instead).
2754 See the related functions putchunk, set and calc.
2755 """
2756 return self._swigobj.putregion(pixels, pixelmask, region, list, usemask, locking, replicate)
2758 def rebin(self, outfile, bin, region='', mask='', dropdeg=False, overwrite=False, stretch=False, crop=False):
2759 """This application rebins the current image by the specified integer binning
2760 factors for each axis. It supports both float valued and complex valued images.
2761 The corresponding output pixel value is the average of the
2762 input pixel values. The output pixel will be masked bad if there
2763 were no good input pixels. A polarization axis cannot be rebinned.
2765 The binning factors array must contain at least one element and no more
2766 elements than the number of input image axes. If the number of elements
2767 specified is less than the number of image axes, then the remaining axes
2768 not specified are not rebinned. All specified values must be positive. A
2769 value of one indicates that no rebinning of the associated axis will occur.
2771 Binning starts from the origin pixel of the bounding box of the selected region or
2772 the origin pixel of the input image if no region is specified. The value of crop
2773 is used to determine how to handle cases where there are pixels
2774 at the end of the axis that do not form a complete bin. If crop=True,
2775 extra pixels at the end of the axis are discarded. If crop=False, the remaining
2776 pixels are averaged into the final bin along that axis. Should the length
2777 of the axis to be rebinned be an integral multiple of the associated binning
2778 factor, the value of crop is irrelevant.
2780 A value of dropdeg=True will result in the output image not containing
2781 axes that are degenerate in the specified region or in the input image if no
2782 region is specified. Note that, however, the binning
2783 factors array must still account for degenerate axes, and the binning
2784 factor associated with a degenerate axis must always be 1.
2786 If {stfaf outfile} is given, the image is written to the specified
2787 disk file. If {stfaf outfile} is unset, the Image tool is
2788 associated with a temporary image. This temporary image may be in
2789 memory or on disk, depending on its size. When you destroy the
2790 on-the-fly Image tool returned by this function (with the done function) this
2791 temporary image is deleted.
2792 """
2793 return _wrap_image(swig_object=self._swigobj.rebin(outfile, bin, region, mask, dropdeg, overwrite, stretch, crop))
2795 def regrid(self, outfile='', shape=[ int(-1) ], csys={ }, axes=[ int(-1) ], region={ }, mask='', method='linear', decimate=int(10), replicate=False, doref=True, dropdeg=False, overwrite=False, force=False, asvelocity=False, stretch=False):
2796 """This function regrids the current image onto a grid specified by the given
2797 coordinate system. The shape of the output image may also be specified.
2799 The coordinate system must be specified via a cs tool (using cs.torecord()). It
2800 is optional; if not specified, the coordinate system from the input image (ie,
2801 the one to which you are applying the regrid function) is used. The order of the
2802 coordinates and axes in the output image is always the same as the input image.
2803 It simply finds the relevant coordinate in the supplied coordinate system in
2804 order to determine the regridding parameters. The supplied coordinate system
2805 must have at least as many coordinates as are required to accomodate the axes
2806 that are to be regridded (eg, if the first two axes are to be regridded, and
2807 these belong to a direction coordinate, one direction coordinate in the supplied
2808 coordinate system is required). Coordinates pertaining to axes that are not
2809 being regridded are supplied from the input image, not the specified coordinate
2810 system.
2812 Reference changes are handled (eg, J2000 to B1950, LSR to TOPO). In general, the
2813 conversion machinery attempts to work out how sophisticated it needs to be (eg,
2814 is the regridding being done from LSR to LSR or from LSR to TOPO). However, it
2815 errs on the side of caution if the conversion machine requires more information
2816 than it actually needs. For full frame conversions, one needs to know things
2817 like the position on the earth's surface where the observation occurred,
2818 direction (celestial coordinates) of observation, and time of observation.
2820 If you get such errors and you are doing a frame conversion, then that means you
2821 must insert some extra information into the coordinate system of your image.
2822 Most likely it's the time (in which case you can use cs.setepoch()) and/or
2823 position (in which case you can use cs.settelescope()) that are missing. If you
2824 get these errors and you are certain that you are not specifying a frame change
2825 (eg, regrid LSR to LSR) then try setting doref=False. This will (silently)
2826 bypass all possible frame conversions. Note that if you are requesting a frame
2827 conversion and you set doref=False, no warnings will be emitted and the output
2828 image will likely be nonsensical.
2830 If you regrid a plane holding a direction coordinate and the units are Jy/pixel,
2831 then the output is scaled to conserve flux (roughly; just one scale factor at
2832 the reference pixel is computed).
2834 Regridding of complex-valued images is supported. The real and imaginary parts
2835 are regridded independently and the resulting regridded pixel values are
2836 combined to form the regridded, complex-valued image.
2838 A variety of interpolation schemes are provided (you need only specify the first
2839 three characters to the method parameter). The cubic interpolation is
2840 substantially slower than linear, and often the improvement is modest. By
2841 default, linear interpolation is used.
2843 You specify the shape of the output image using the shape parameter and which
2844 output axes you want to regrid. Note that a Stokes axis cannot be regridded
2845 (you will get a warning if you try).
2847 The axes parameter cannot be used to discard axes from the output image; it can
2848 only be used to specify which output axes are going to be regridded and which
2849 are not. Any axis that you are not regridding must have the same output shape as
2850 the input image shape for that axis.
2852 The axes parameter can also be used to specify the order in which the output
2853 axes are regridded. This may give you significant performance benefits. For
2854 example, imagine we are going to regrid a spectral-line cube of shape
2855 [512,512,1204] to shape [256,256,32]. If you specified axes=[0,1,2] then first,
2856 the direction axes would be regridded for each of the 1024 pixels (and stored in
2857 a temporary image). Then each spectral profile at each spatial location in the
2858 temporary image would be regridded to 32 pixels. You could speed this process
2859 up significantly by setting axes=[2,0,1]. In this case, first each spectral
2860 profile would be regridded to 32 pixels, and then each plane of the 32 pixels
2861 would be regridded. Note that the order of axes does not affect the order of the
2862 shape parameter; ie, it should be given in the natural pixel axis order of the
2863 image ()[256,256,32] in both cases of this example).
2865 You can also specify a region to be applied to the input image. If you do this,
2866 you need to be careful with the output shape for non-regridded axes (must match
2867 that of the region - use function ia.boundingbox() to determine that).
2869 If the outfile parameter is specified, the image is written to the specified
2870 disk file. If this parameter is unset, the on-the-fly image analysis tool
2871 returned by this method is associated with a temporary image. This temporary
2872 image may be in memory or on disk, depending on its size. When you destroy the
2873 on-the-fly image analysis tool (with either the ia.close() or ia.done()
2874 methods), this temporary image is deleted.
2876 The replicate parameter can be used to simply replicate pixels rather than
2877 regridding them. Normally replicate=False, for every output pixel, its world
2878 coordinate is computed and the corresponding input pixel found (then a little
2879 interpolation grid is generated). If replicate=True, then for every output
2880 axis, a vector of regularly sampled input pixels is generated (based on the
2881 ratio of the output and input axis shapes). So this just means the pixels get
2882 replicated (by whatever interpolation scheme you use) rather than regridded in
2883 world coordinate space. This process is much faster, but its not a true world
2884 coordinate based regrid.
2886 As decribed above, when replicate=False, a coordinate is computed for each
2887 output pixel; this is an expensive operation. The decimate parameter allows you
2888 to decimate the computation of that coordinate grid to a sparse grid, which is
2889 then filled in via fast interpolation. The default is decimate=10. The number
2890 of pixels per axis in the sparse grid is the number of output pixels for that
2891 axis divided by the decimation factor. A factor of 10 does pretty well. You may
2892 find that for very non-linear coordinate systems (e.g. very close to the pole)
2893 that you have to reduce the decimation factor. You may also have to reduce the
2894 decimation factor if the number of pixels in the output image along an axis to
2895 be regridded is less than about 50, or the output image may be completely
2896 masked.
2898 If one of the axes to be regridded is a spectral axis and asvelocity=True, the
2899 axis will be regridded to match the velocity, not the frequency, description of
2900 the template spectral coordinate. Thus the output pixel values will correspond
2901 only to the velocity, not the frequency, of the output axis.
2903 Sometimes it is useful to drop axes of length one (degenerate axes). Setting
2904 the dropdeg parameter to True will do that. It will discard the axes from the
2905 input image. Therefore the output shape and coordinate system that you supply
2906 must be consistent with the input image after the degenerate axes are dropped.
2908 The force parameter can be used to force all specified axes to be regridded,
2909 even if the algorithm determines that they don't need to be (because the input
2910 and output coordinate information is identical).
2912 The cs tool has a useful method, cs.setreferencelocation(), that can be used to
2913 keep a specific world coordinate in the center of an image when regridding (see
2914 example below).
2916 The output pixel mask will be True (good) unless the regridding failed to find a
2917 value for that output pixel in which case it will be False. For example, if the
2918 total input mask (default input pixel mask plus OTF mask) for all of the
2919 relevant input pixels were masked bad then the output pixel would be masked
2920 (False).
2922 MULTIPLE AXIS COORDINATES LIMITATION
2923 Some cooordinates contain multiple axes. For example, a direction coordinate
2924 contains both longitude-like and latitude-like axes. A linear coordinate can
2925 also contain multiple axes. When you regrid *any* axis from a coordinate which
2926 contains multiple axes, you must fully specify the coordinate information for
2927 all axes in that coordinate in the coordinate system that you provide. For
2928 example, if you have a linear coordinate with two axes and you want to regrid
2929 axis one only. In the coordinate system you provide, the coordinate information
2930 for axis two (not being regridded) must correctly be a copy from the input
2931 coordinate system (it won't be filled in for you).
2933 If an image has per-plane beams and one attempts to regrid the spectral axis,
2934 an exception is thrown.
2936 IMPORTANT NOTE ABOUT FLUX CONSERVATION
2937 in general regridding is inaccurate for images that the angular resolution is
2938 poorly sampled. A check is done for such cases and a warning message is emitted
2939 if a beam present. However, no such check is done if there is no beam present.
2940 To add a restoring beam to an image, use ia.setrestoringbeam().
2941 """
2942 return _wrap_image(swig_object=self._swigobj.regrid(outfile, shape, csys, axes, region, mask, method, decimate, replicate, doref, dropdeg, overwrite, force, asvelocity, stretch))
2944 def transpose(self, outfile='', order=[ ]):
2945 """This method transposes the axes in the input image to the specified
2946 order. The associated pixel and mask values and coordinate system are transposed.
2948 If the outfile parameter is empty, only a temporary image is created; no output image
2949 is written to disk.
2951 The order parameter describes the mapping of the input axes to the output axes.
2952 It can be one of three types: a non-negative integer, a string, or a list of
2953 strings. If a string or non-negative integer, it should contain
2954 zero-based digits describing the new order of the input axes. It must
2955 contain the same number of (unique) digits as the number of input axes. For example,
2956 specifying reorder="1032" or reorder=1032 for a four axes image maps input axes
2957 1, 0, 3, 2 to output axes 0, 1, 2, 3. In the case of order being a nonnegative integer
2958 and the zeroth axis in the input being mapped to zeroth axis in the output, the zeroth
2959 digit is implicitly understood to be 0 so that to transpose an image where one would
2960 use a string order="0321", one could equivalently specify an int order=321.
2961 IMPORTANT: When specifying a non-negative integer and mapping the zeroth axis of
2962 the input to the zeroth axis of the output, do *not* explicitly specify the leading
2963 0; eg, specify order=321 rather than order=0321. Python interprets an integer with
2964 a leading 0 as an octal number.
2966 Because of ambiguity for axes numbers greater than nine, using string or integer order
2967 specifications cannot handle images containing more than 10 axes.
2968 The order parameter can also be specified as a list of strings which uniquely minimally match,
2969 ignoring case, the image axis names (ia.coordsys().names()).
2970 So to reorder an image with right ascension, declination, and frequency axes, one could
2971 specify order=["d", "f", "r"] or equivalently ["decl", "frequ", "right a"]. Note that
2972 specifying "ra" for the right ascension axis will result in an error because "ra" does
2973 not match the first two characters of right ascension.
2974 Axes can be simultaneously inverted in cases where order is a string or an array of
2975 strings by specifying negative signs in front of the axis/axes to be inverted. So,
2976 in a 4-D image, order="-10-3-2" maps input axes 1, 0, 3, 2 to output axes 0, 1, 2, 3
2977 and reverses the direction and values of input axes 1, 3, and 2.
2979 """
2980 return _wrap_image(swig_object=self._swigobj.transpose(outfile, order))
2982 def rotate(self, outfile='', shape=[ int(-1) ], pa=[ ], region={ }, mask='', method='cubic', decimate=int(0), replicate=False, dropdeg=False, overwrite=False, stretch=False):
2983 """This function rotates two axes of an image. These axes are either
2984 those associated with a Direction coordinate or with a Linear
2985 coordinate. The Direction coordinate takes precedence.
2986 If rotating a Linear coordinate, it must hold precisely two axes.
2988 The method is that the Coordinate is rotated and then the input
2989 image is regridded to the rotated Coordinate System.
2991 If the image brightness units are Jy/pixel then the output is scaled to
2992 conserve flux (roughly; just one scale factor at the reference pixel is
2993 computed).
2995 A variety of interpolation schemes are provided (you need only specify
2996 the first three characters to {stfaf method}). The cubic
2997 interpolation is substantially slower than linear. By default you get
2998 cubic interpolation.
3000 You can specify the shape of the output image ({stfaf shape}).
3001 However, all axis that are not regrided retain the same output shape
3002 as the input image shape for that axis. Only the direction coordinate
3003 axes are regridded.
3005 You can also specify a region to be applied to the input image. If
3006 you do this, you need to be careful with the output shape for
3007 non-regridded axes (must match that of the region - use function
3008 boundingbox to find that out).
3010 If {stfaf outfile} is given, the image is written to the specified
3011 disk file. If {stfaf outfile} is unset, the on-the-fly Image tool
3012 returned by this function is associated with a temporary image. This
3013 temporary image may be in memory or on disk, depending on its size.
3014 When you destroy the on-the-fly Image tool (with the done function) this
3015 temporary image is deleted.
3017 The argument {stfaf replicate} can be used to simply replicate pixels
3018 rather than regridding them. Normally ({stfaf replicate=F}), for every
3019 output pixel, its world coordinate is computed and the corresponding
3020 input pixel found (then a little interpolation grid is generated). If
3021 you set {stfaf replicate=T}, then what happens is that for every output
3022 axis, a vector of regularly sampled input pixels is generated (based on
3023 the ratio of the output and input axis shapes). So this just means the
3024 pixels get replicated (by whatever interpolation scheme you use) rather
3025 than regridded in world coordinate space. This process is much faster,
3026 but its not a true world coordinate based regrid.
3028 As decribed above, when {stfaf replicate} is False, a coordinate is
3029 computed for each output pixel; this is an expensive operation. The
3030 argument {stfaf decimate} allows you to decimate the computation of
3031 that coordinate grid to a sparse grid, which is then filled in via
3032 fast interpolation. The default for {stfaf decimate} is 0 (no
3033 decimation). The number of pixels per axis in the sparse grid is the
3034 number of output pixels for that axis divided by the decimation
3035 factor. A factor of 10 does pretty well. You may find that for very
3036 non-linear coordinate systems (e.g. very close to the pole) that you
3037 have to reduce the decimation factor.
3039 The output pixelmask will be good (T) unless the regridding failed to
3040 find a value for that output pixel in which case it will be bad (F).
3041 For example, if the total input mask (default input pixelmask plus OTF
3042 mask) for all of the relevant input pixels were masked bad
3043 then the output pixel would be masked bad (F).
3044 """
3045 return _wrap_image(swig_object=self._swigobj.rotate(outfile, shape, pa, region, mask, method, decimate, replicate, dropdeg, overwrite, stretch))
3047 def rotatebeam(self, angle=[ ]):
3048 """This method rotates the attached image's beam(s) counterclockwise through the specified angle.
3049 This is the same thing as increasing the position angle(s) of the beam(s) by the specified angle.
3050 If the image does not have a beam, no changes to the image are made. If the image has multiple
3051 beams, all the beams are rotated through the same angle.
3052 """
3053 return self._swigobj.rotatebeam(angle)
3055 def rename(self, name, overwrite=False):
3056 """This function renames the imagefile associated with the imagetool.
3057 If a file with name {stfaf name} already exists, you can overwrite it
3058 with the argument {stfaf overwrite}; otherwise a fail will
3059 result.
3060 """
3061 return self._swigobj.rename(name, overwrite)
3063 def replacemaskedpixels(self, pixels=[ ], region={ }, mask='', update=False, list=False, stretch=False):
3064 """This application replaces the values of all pixels whose total input mask
3065 (default input pixelmask and OTF mask) is bad (F) with the specified
3066 value. It supports both float valued and compplex valued images.
3068 If the argument {stfaf update} is F (the default), the actual pixelmask
3069 is left unchanged. That is, masked pixels remain masked. However, if
3070 you set {stfaf update=T} then the pixelmask will be updated so that the
3071 pixelmask will now be T (good) where the {bf total} input mask was F
3072 (bad).
3074 See maskhandler for information
3075 on how to set the default pixelmask.
3077 There are a few ways in which you can specify what to replace the
3078 masked pixel values by.
3080 begin{itemize}
3082 item First, you can give the {stfaf pixels} argument a simple numeric
3083 scalar (e.g. {cf pixels=1.0}). Then, all masked values will be
3084 replaced by that value.
3086 item Second, you can give a scalar
3087 htmladdnormallink{LEL}{../../notes/223/223.html} expression string
3088 (e.g. {cf pixels='min(myimage)'}). Then, all masked values will be
3089 replaced by the scalar that results from the expression. If the scalar expression
3090 is illegal (e.g. in the expression {cf pixels='min(myimage)'} there
3091 were no good pixels in {sff myimage}) then the value 0 is used for
3092 replacement.
3094 item Third, you can give a
3095 htmladdnormallink{LEL}{../../notes/223/223.html} expression string
3096 which has the same shape as the imagefile you are applying the
3097 function to. For example, putting {cf pixels='myotherimage'} means
3098 replace all masked pixels in this imagefile with the equivalent pixel
3099 in the imagefile called {sff myotherimage}.
3101 Your expression might be quite complex, and you can think of it as
3102 producing another masked lattice. However, in the replace process, the
3103 mask of that expression lattice is ignored. Thus, only the mask of
3104 the imagefile you are replacing and the pixel values of the expression
3105 lattice are relevant.
3107 The expression must conform with the subimage formed by applying the
3108 region to the image (i.e. that associated with this Image tool). If
3109 you use the {stfaf mask} argument as well, the region is applied to
3110 it as well (see examples).
3112 end{itemize}
3113 """
3114 return self._swigobj.replacemaskedpixels(pixels, region, mask, update, list, stretch)
3116 def beamarea(self, channel=int(-1), polarization=int(-1)):
3117 """Get the area of the image's restoring beam. If a non-negative channel and non-negative polarization
3118 are specified, the area for the beam associated with that channel and polarization will be
3119 returned. The return value will be a dictionary containing the keys 'arcsec2' and 'pixels', and
3120 the associated values will be the beam area in arcsec2 and in pixels, respectively. If both
3121 channel and polarization are set to negative values, then a dictionary with the same keys will
3122 be returned, and the values will be either scalars (if the image has a single traditional
3123 beam) or arrays if the image has multiple beams. In the latter case, the arrays will have
3124 shapes indicative of the number of channels and number of polarizations. If the image has
3125 a spectral axis but not a polarization axis, the returned arrays will be one dimensional and
3126 have a length equal to the number of channels. Similarly, if the image has a polarization
3127 axis but not a spectral axis, the arrays will be one dimensional and have a lenghts equal
3128 to the number of polarizations. If the image has both a spectral and polarization axis,
3129 the returned arrays will be two dimensional with shape (m, n) where m is the length
3130 of the first channel or polarization axis, and n is the length of the second channel or
3131 polarization axis. So, if an image has shape [200, 200, 10, 4] with 10 channels and
3132 4 stokes, the returned arrays will have shapes of (10, 4) representing the spectral
3133 axis as the first axis and the polarization axis as the second. If, instead, the image
3134 shape is [200, 200, 4, 10] again with 10 channels and 4 stokes, the shape of the
3135 returned arrays will be (4, 10) since the polarization axis precedes the spectral
3136 axis.
3138 """
3139 return self._swigobj.beamarea(channel, polarization)
3141 def restoringbeam(self, channel=int(-1), polarization=int(-1), mbret='list'):
3142 """This function returns the restoring beam(s), if any, of the attached image.
3144 If mbret="matrix" (or "m"), an exception will be thrown if the attached image
3145 does not have per-plane beams, or if channel or polarization is specified to be
3146 non-negative. See below for more information on the mbret parameter.
3148 If mbret="list" (or "l") and the attached image has no restoring beam(s), an
3149 empty dictionary is returned.
3151 If the image has a single traditional restoring beam and mbret="list",
3152 the beam is returned as a dictionary no matter what channel and polarization
3153 values are specified. This dictionary has keys "major",
3154 "minor", andi "positionangle", and each of these fields
3155 contains a dictionary with "value" and "unit" keys which
3156 provide the quantity associated with each beam parameter.
3158 If mbret="list", and if the image has per-plane beams, and the image
3159 has both spectral and polarization axes, and if both channel and polarization
3160 are set to non-negative values that are less than the number of planes along
3161 each parameter's representative axis, the beam for that particular
3162 channel/polarization pair is returned. If at least one of these values is set
3163 to greater or equal to the number of planes along that parameter's
3164 representative axis, and that axis exists and is not degenerate, then an
3165 exception will be thrown. In the case of a non-extant axis or a degenerate axis,
3166 the parameter associated with that axis can be set to any value and a beam will
3167 be returned, because in that case there is exactly one beam for each plane
3168 along the other parameter's representative axis, and so there is no ambiguity
3169 regarding which beam to return. In the case where both spectral and polarization
3170 axes exist and are not degenerate, an exception will be thrown if one of channel
3171 or polarization is set to a non-negative value and the other is set to a
3172 negative value. In the above cases in which a beam is returned, the returned
3173 dictionary will have the same structure as that previously described in the
3174 single beam case.
3176 If the image contains multiple beams and both channel and polarization are
3177 specified to be negative, the structure of the returned dictionary depends on
3178 the specified value of mbret. Supported values of this parameter are
3179 "list" and "matrix" (case insensitive, mimimum match
3180 supported). In both cases, the returned dictionary will contain the keys
3181 "nChannels", which contains an integer value equal to the number of
3182 channels, and "nStokes", which contains an integer value equal to the
3183 number of polarizations. In the case where one axis doesn't exist, the
3184 associated value will be 1.
3186 In the case of mbret="list", the returned dictionary will contain the
3187 key "beams", which contains a sub-dictionary of information for all
3188 beams. This subdictionary contains keys "*0" through
3189 "*(c-1)", where c is the value associated with the
3190 "nChannels" key. Each of these keys has an associated value which is a
3191 subdictionary containing keys "*0" through "*(p-1)", where p
3192 is the value associated associated with the "nStokes" key. Each
3193 of these keys has an assciated value of a beam dictionary, the structure of
3194 which is the same as described previously for a single beam image, which is
3195 associated with that particular channel/polarization pair.
3197 In the case of mbret="matrix", the returned dictionary will have keys
3198 "major", "minor", and "pa" which hold values
3199 representing the major axes, the minor axes, and the position angle,
3200 respectively, of all beams. Each of these fields contains a dictionary which
3201 represents a quantity with keys "unit" and "value". The
3202 "unit" field contains a string representing the unit of the
3203 associated value. The "value" field contains a matrix with values
3204 corresponding to the associated parameter for the associated beam. In the case
3205 where an image does not have both a spectral and polarization axis, these
3206 matrices will have shape (n, 1), where n is the number of planes on the extant
3207 axis. In the case where an image has both spectral and polarization axes, the
3208 matrices will have shape (m, n) where m is the number of channels if the
3209 spectral axis precedes the polarization axis, or the number of polarization
3210 planes otherwise, and n is the number of planes on the axis not represented by
3211 m. So, an image with a spectral axis that has 5 channels and precedes the
3212 polarization axis that has 4 planes, the returned matrix shape will be (5, 4).
3213 In the case where the number of spectral and polarization planes are the same
3214 as in the previous example, but the polarization axis precedes the spectral
3215 axis, the matrix shape will be (4, 5). Thus, the matrix shape pair follows the
3216 order of the image spectral and polarization axes. This mode is useful for
3217 returning to the python UI the major axes, minor axes, and position angles
3218 directly as numpy arrays.
3220 The restoring beam(s) of an image may be set with image tool method
3221 setrestoringbeam.
3222 """
3223 return self._swigobj.restoringbeam(channel, polarization, mbret)
3225 def sepconvolve(self, outfile='', axes=[ int(-1) ], types=[ '' ], widths=[ ], scale=float(-1), region={ }, mask='', overwrite=False, stretch=False):
3226 """This function does Fourier-based
3227 convolution of an imagefile by a specified separable kernel.
3229 If {stfaf outfile} is given, the image is written to the specified
3230 disk file. If {stfaf outfile} is unset, the on-the-fly Image tool
3231 returned by this function is associated with a temporary image. This
3232 temporary image may be in memory or on disk, depending on its size.
3233 When you destroy the Image tool (with the done function) this
3234 temporary image is deleted.
3236 You specify which axes of the image you wish to convolve, by what kernel
3237 of what width. The kernel types can be shortened to {cf 'gauss',
3238 'hann'} and {cf 'box'}.
3240 You specify the widths of the convolution kernels via the argument
3241 {stfaf widths}. The values can be specified as a vector of three
3242 different types.
3244 begin{itemize}
3246 item Quantity - for example {stfaf widths=qa.quantity("1arcsec 0.00001rad")}.
3247 Note that you can use pixel units, viz. {stfaf widths=qa.quantity("10pix 0.00001rad")}
3248 see below.
3250 item String - for example {stfaf widths="1km 2arcsec"} (i.e. a string that
3251 qa.quantity() accepts).
3255 item Numeric - for example {stfaf widths=[10,20]}. In this case,
3256 the units of the widths are assumed to be in pixels.
3258 end{itemize}
3260 The interpretation of {stfaf widths} depends upon the kernel type.
3262 begin{itemize}
3264 item Gaussian - the specified width is the full-width at
3265 half-maximum.
3267 item Boxcar (tophat) - the specified width is
3268 the full width.
3270 item Hanning - The kernel is $z[i] = 0.25*y[i-1] + 0.5*y[i] +
3271 0.25*y[i+1]$. The width is always 3 pixels, regardless of what
3272 you give (but you still have to give it !).
3274 end{itemize}
3276 The scaling of the output image is determined by the argument {stfaf scale}.
3277 If you leave it unset, then autoscaling will be invoked which means that
3278 the convolution kernels will all be normalized to have unit volume
3279 to as to conserve flux.
3281 If you do not leave {stfaf scale} unset, then the convolution kernel
3282 will be scaled by this value (it has peak unity before application
3283 of this scale factor).
3285 Masked pixels will be assigned the value 0.0 before convolution.
3286 The output mask is the combination (logical OR) of the default input
3287 pixelmask (if any) and the OTF mask. Any other input pixelmasks
3288 will not be copied. Use function
3289 maskhandler if you need to copy other
3290 masks too.
3292 See also the other convolution functions
3293 convolve2d,
3294 convolve and
3295 hanning.
3296 """
3297 return _wrap_image(swig_object=self._swigobj.sepconvolve(outfile, axes, types, widths, scale, region, mask, overwrite, stretch))
3299 def set(self, pixels=[ ], pixelmask=int(-1), region={ }, list=False):
3300 """This function replaces data and/or mask values within the image in the
3301 specified region. You can think of it as a simplified version of the
3302 image calculator.
3304 Unlike the calc function, you can
3305 only set a scalar value for all pixels in the specified region. For
3306 example, it can be useful to set a whole image to one value, or a mask
3307 in a region to one value. Although you could do that with the related
3308 functions putregion and
3309 putchunk, you would have to make an
3310 array of the shape of the image and if that is large, it could be
3311 resource expensive.
3313 The value for the pixels is specified with the {stfaf pixels} argument. It can
3314 be given as either a Lattice Expression Language (or LEL) expression
3315 string or a simple numeric scalar. See htmladdnormallink{note
3316 223}{../../notes/223/223.html} for a detailed description of the LEL
3317 expression syntax. If you give a LEL expression it must be a scalar
3318 expression.
3320 Note that any default mask is {em ignored} by this function when you
3321 set pixel values. This is different from
3322 calc where the extant mask is
3323 honoured.
3325 The value for the pixel mask is specified with the {stfaf pixelmask}
3326 argument ({cf T, F, unset}). If it's {cf unset} then the mask is not
3327 changed.
3329 If you specify {stfaf pixelmask=} T or F, then the mask that is affected is
3330 the current default mask (see
3331 maskhandler). If there is no mask, a
3332 mask is created for you and made the default mask.
3333 """
3334 return self._swigobj.set(pixels, pixelmask, region, list)
3336 def setbrightnessunit(self, unit):
3337 """This function sets the image brightness unit. Both float and complex
3338 valued images are supported.
3339 You can get the brightness unit with function
3340 brightnessunit.
3341 """
3342 return self._swigobj.setbrightnessunit(unit)
3344 def setcoordsys(self, csys):
3345 """This method replaces the coordinate system in the image. Coordinate systems are
3346 manipulated with the cs (coordintate system) tool. The coordinate system can be
3347 recovered from an image via the coordsys() method of the image analysis (ia)
3348 tool.
3350 Note that changing the coordinate system using the cs tool has no effect on the
3351 original image, until it is replaced with this method; the value returned
3352 by coordsys() is a copy of, not a reference to, the image's coordinate system.
3353 """
3354 return self._swigobj.setcoordsys(csys)
3356 def sethistory(self, origin='', history=[ ]):
3357 """A casa imagefile can accumulate history information
3358 from an input fits file or by you writing something into
3359 it explicitly with this function. Each element of
3360 the input vector is one line of history. The new
3361 history is appended to the old.
3364 You can recover the history information with function
3365 history.
3366 """
3367 return self._swigobj.sethistory(origin, history)
3369 def setmiscinfo(self, info):
3370 """A CASA image can include user-specified or miscellaneous metadata. This metadata
3371 is stored in a data structure referred to as a miscinfo record. For example, the
3372 FITS reader ia.fromfits() puts header keywords it doesn't otherwise use into
3373 such a record. The miscinfo record is not required to be populated, though.
3375 This method sets the miscinfo record of an image. Note that this method will
3376 overwrite, not add to, an existing miscinfo record. Thus if the user wishes
3377 to augment an existing record, they must first capture the existing record using
3378 the ia.miscinfo() method, modify the captured record, and then replace the
3379 record in the image using setmiscinfo() by passing it the modified record. The
3380 FITS writer will attempt to write all the fields in the miscinfo record to the
3381 FITS file it creates. It can do so for scalars and 1-dimensional arrays. Records
3382 will be omitted, and multi-dimensional arrays will be flattened into
3383 one dimensional arrays.
3384 """
3385 return self._swigobj.setmiscinfo(info)
3387 def shape(self):
3388 """The shape of an image is a vector holding the length of each axis of
3389 the image. Although this information is also available in the summary function, it is
3390 so useful that it can be obtained directly. Both Float and Complex valued
3391 images are supported.
3392 """
3393 return self._swigobj.shape()
3395 def setrestoringbeam(self, major=[ ], minor=[ ], pa='Not specified', beam={ }, remove=False, log=True, channel=int(-1), polarization=int(-1), imagename=''):
3396 """This method sets the restoring beam(s) for an image or removes an existing beam(s).
3398 An image must have exactly one of the following states:
3400 1. An image can have a single "traditional" or global beam. In that case, the
3401 beam applies to every channel and polarization in the image.
3403 2. If an image has more than one spectral channel or more than one polarization,
3404 it can have a set of beams. In this case, each channel and/or polarization will
3405 have its own beam.
3407 3. An image can have neither a global beam nor a beam set.
3409 It is never permissible for an image to have both a traditional (global) beam
3410 and a set of per-plane beams. Task and method behavior is undefined in that case
3411 and any resulting products are considered corrupt.
3413 RULES FOR BEAM MODIFICATION
3415 1. If remove=true, any existing beam(s) are removed.
3417 2. Else if imagename is specified, the beam(s) from the specified image will
3418 be copied to the image being accessed by the image tool. Multiple beams
3419 may be copied, but the two images must have the same number of frequency
3420 channels and polarization planes. If not, an exception is thrown.
3422 3. Else if the beam parameter is specified, it will be used. It must be fully
3423 specified. It must have exactly three items with keys "major",
3424 "minor", and "positionangle". Each of these keys
3425 must be a proper quantity dictionary with keys "value" and
3426 "unit". The units for all three items should be angular. The
3427 "major" and "minor" items must have non-negative values.
3428 If any of these conditions is not met, an exception will be thrown.
3430 4. Else the "major", "minor", and "pa"
3431 parameters must be specified. Any or all of these may be a quantity string
3432 (eg, "2arcsec", a quantity dictionary with keys "value"
3433 and "unit". or a numerical value. "major" and
3434 "minor" must have non-negative values. If any of these conditions
3435 is not met, an exception will be thrown. In the case of a quantity
3436 being specified the "unit" should be an angular unit or no unit.
3437 In the case of no unit or a numerical value specified, if the image already
3438 has a beam then the corresponding unit for that parameter in the current beam
3439 will be used. If the image has no beam, then "arcsec" is used for
3440 "major" and "minor", and "deg" is used for
3441 "pa".
3443 If an image has no beams, a traditional (global) beam can be added by setting
3444 both channel and polarization to negative values.
3446 If an image has no beams, a set of per-plane beams can be added by setting
3447 either or both channel and/or polarization to a non-negative value. In this
3448 case, a number of per-plane beams are added consistent with the image and they
3449 are all set to be the same with parameters equal to those specified by either
3450 the beam or major/minor/pa parameters.
3452 If an image has a traditional beam, it can be modified by setting both channel
3453 and polarization to negative values. If one or both is not set to a negative
3454 value, an exception is thrown, and nothing is modified.
3456 If an image has a set of per plane beams, one at a time of these can be modified
3457 by setting the appropriate channel number and/or polarization number. All the
3458 per-plane beams can be modified to the same values in one go by setting both
3459 channel and polarization to negative values. Also, in the case where an image
3460 has multiple channels, the beams associated with all channels for a given
3461 polarization can be modified to the same beam by setting polarization equal to
3462 the desired polarization plane number and by setting channel to a negative
3463 value. Similarly, in the case where an image has multiple polarizations, the
3464 beams associated with all polarizations for a given spectral channel can be
3465 modified to the same beam by setting channel equal to the desired spectral
3466 channel number and by setting polarization to a negative value.
3467 """
3468 return self._swigobj.setrestoringbeam(major, minor, pa, beam, remove, log, channel, polarization, imagename)
3470 def statistics(self, axes=[ int(-1) ], region='', mask='', includepix=[ float(-1) ], excludepix=[ float(-1) ], list=False, force=False, disk=False, robust=False, verbose=False, stretch=False, logfile='', append=True, algorithm='classic', fence=float(-1), center='mean', lside=True, zscore=float(-1), maxiter=int(-1), clmethod='auto', niter=int(3)):
3471 """This method computes statistics from the pixel values in the image. You can
3472 then list them and retrieve them (into a record) for further analysis. This
3473 method supports only real valued images.
3475 The names of the fields in the returned record are summarized below:
3477 - npts: the number of unmasked points used
3479 - sum: the sum of the pixel values: $sum I_i$
3481 - flux: flux or flux density, see below for details
3483 - sumsq: the sum of the squares of the pixel values: $sum I_i^2$
3485 - mean: the mean of pixel values: $bar{I} = sum I_i / n$
3487 - sigma: the standard deviation about the
3488 mean: $sigma^2 = (sum I_i - bar{I})^2 / (n-1)$
3490 - rms: the root mean square: $sqrt {sum I_i^2 / n}$
3492 - min: minimum pixel value
3494 - max: the maximum pixel value
3496 - median: the median pixel value (if {stfaf robust=T})
3498 - medabsdevmed: the median of the absolute deviations from the median
3499 (if {stfaf robust=T})
3501 - quartile: the inter-quartile range (if {stfaf
3502 robust=T}). Find the points which are 25% largest and 75% largest
3503 (the median is 50% largest).
3505 - q1: The first quartile. Reported only if robust=T.
3507 - q3: The third quartile. Reported only if robust=T.
3509 - blc: the absolute pixel coordinate of the bottom left
3510 corner of the bounding box of the region of interest. If 'region' is
3511 unset, this will be the bottom left corner of the whole image.
3513 - blcf: the formatted absolute world coordinate of the bottom left corner of the bounding box of the region of interest.
3515 - trc: the absolute pixel coordinate of the top right corner of the bounding box of the region of interest.
3517 - trcf: the formatted absolute world coordinate of the top right corner of the bounding box of the region of interest.
3519 - minpos: absolute pixel coordinate of minimum pixel value
3521 - maxpos: absolute pixel coordinate of maximum pixel value
3523 - minposf: formatted string of the world coordinate of
3524 the minimum pixel value
3526 - maxposf: formatted string of the world coordinate of
3527 the maximum pixel value
3530 The last four fields only appear if you evaluate the statistics over all
3531 axes in the image. As an example, if the returned record is captured in
3532 {stfaf 'mystats'}, then you could access the 'mean' field via
3533 {cf print mystats['mean']}.
3535 If there are no good points (e.g. all pixels are masked bad in the
3536 region), then the length of these fields will be 0 (e.g. {cf
3537 len(mystats['npts'])==0}).
3539 You have no control over which statistics are listed to the logger,
3540 you always get the same selection. You can choose to list the
3541 statistics or not (argument {stfaf list}).
3543 As well as the simple (and faster to calculate) statistics like means
3544 and sums, you can also compute some robust (quantile-like) statistics. Currently
3545 these are the median, median absolute deviations from the median,
3546 the first and third quartiles, and the inner-quartile range. Because these are computationally
3547 expensive, they are only computed if robust=True.
3549 Note that if the axes are set to all of the axes in the image (which is
3550 the default) there is just one value per statistic.
3552 You have control over which pixels are included in the statistics computations
3553 via the {stfaf includepix} and {stfaf excludepix} arguments. These vectors
3554 specify a range of pixel values for which pixels are either included or
3555 excluded. They are mutually exclusive; you can specify one or the
3556 other, but not both. If you only give one value for either of these,
3557 say {stfaf includepix=b}, then this is interpreted as {stfaf
3558 includepix=[-abs(b),abs(b)]}.
3560 This function generates a 'storage' lattice, into which the statistics
3561 are written as a means to improve performance on successive identical runs.
3562 After the initial execution of this method, it is only regenerated as
3563 necessary. For example, if you run the method twice with
3564 identical arguments, the statistics will be directly retrieved from the
3565 storage lattice the second time, and the statistics will not be recomputed.
3566 However, you can force regeneration of the statistics if you set force=True.
3567 VERY IMPORTANT NOTE. If you have an open image tool on which you've run
3568 statistics(), and
3569 change the pixel values of the opened image via another tool or a task, the
3570 opened image tool has no knowledge that pixel values have been changed, and
3571 so if you run ia.statistics() again with that tool, you will very likely
3572 get incorrect results since they come from the statistics stored from a previous
3573 run. First of all, it is highly discouraged that you do anything like this.
3574 Tools maintain state that only they know about internally; if you change
3575 an image that is already opened with another image tool or task, the general
3576 outcome will be undefined when you run methods on the already opened tool.
3577 However, if for some reason you must do this,
3578 first consider changing your algorithm because there should be no reason
3579 to have to do this. But, if you decide to continue down that path that will
3580 likely lead you over the edge of a cliff, then you can set force=True to
3581 ensure statistics are always recomputed.
3583 The storage medium is either in memory or on disk,
3584 depending upon its size. You can force it to disk if you set {stfaf
3585 disk=T}, otherwise it decides for itself.
3587 CURSOR AXES
3588 The axes parameter allows one to set the cursor axes over which statistics
3589 are computed. For example, consider a 3-dimensional image for which axes=[0,2].
3590 The statistics would be computed for each XZ (axes 0 and 2) plane in the
3591 image. One could then examine those statistics as a function of the Y
3592 (axis 1) axis.
3594 Each statistic is stored in an array in its own field in the returned dictionary.
3595 The dimensionality of these arrays is equal to the number of axes over which the
3596 statistics were not evaluated (called the display axes). For example, if the input
3597 image has four axes, and axes=[0], the output statistic arrays will have three dimensions.
3598 If axes=[0, 1], the output statistic arrays will have two dimensions.
3600 The shape of the output arrays when axes has a positive number of elements is based on
3601 the region selection. If there is no region selection, the shape of the statistic arrays
3602 is just the shape of the image along the display (non-cursor) axes. For example, if the
3603 input image has dimensions of 300x400x4x80 and axes=[0, 1], in the absence of a region
3604 selection, the shape of the output statistic arrays will be 4x80. If there is a region
3605 selection, the shape of the output statistic arrays will be determined by the number of
3606 planes along the display axes chosen in the region selection. For example, continuing with
3607 our example, if axes=[0,1] and region=rg.box([0, 0, 1, 20], [299,399, 2, 60]), the output
3608 statistic arrays will have shapes of 2x41. Only the selected planes will be displayed in the
3609 logger output if verbose=True.
3611 In the case where the image has a pixel mask, and/or the mask parameter is specified, and because
3612 of this specification a plane is entirely masked, this element is included in the statistic arrays
3613 (usually with a value of 0). It is not included in the logger output if verbose=True. One can
3614 exclude such elements from computations on the output arrays by using the numpy.extract() method.
3615 For example, to compute the minimum rms value, not including any fully masked planes, one could
3616 use
3618 stats = ia.statistics(...)
3619 rmsmin = numpy.min(numpy.extract(stats['npts']>0, stats['rms']))
3621 Thus in the computation of rmsmin, only the rms elements are considered which have
3622 associated values of npts that are not zero.
3624 ALGORITHMS
3626 Several types of statistical algorithms are supported:
3628 - classic: This is the familiar algorithm, in which all unmasked pixels, subject to any
3629 specified pixel ranges, are used. One may choose one of two methods, which vary only by
3630 performance, for computing classic statistics, via the clmethod parameter. The "tiled"
3631 method is the old method and is fastest in cases where there are a large number of
3632 individual sets of statistics to be computed and a small number of data points per set.
3633 This can occur when one sets the axes parameter, which causes several individual sets of
3634 statistics to be computed. The "framework" method uses the new statistics framework to
3635 compute statistics. This method is fastest in the regime where one has a small number of
3636 individual sets of statistics to calculate, and each set has a large number of points.
3637 For example, this method is fastest when computing statistics over an entire image in one
3638 go (no axes specified). A third option, "auto", chooses which method to use by predicting
3639 which be faster based on the number of pixels in the image and the choice of the axes
3640 parameter.
3642 - fit-half: This algorithm calculates statistics on a dataset created from real and virtual pixel values.
3643 The real values are determined by the input parameters center and lside. The parameter center
3644 tells the algorithm where the center value of the combined real+virtual dataset should be. Options
3645 are the mean or the median of the input image's pixel values, or at zero. The lside parameter tells
3646 the algorithm on which side of this center the real pixel values are located. True indicates that
3647 the real pixel values to be used are <= center. False indicates the real pixel values to be used
3648 are >= center. The virtual part of the dataset is then created by reflecting all the real values
3649 through the center value, to create a perfectly symmetric dataset composed of a real and a virtual
3650 component. Statistics are then calculated on this resultant dataset. These two parameters are
3651 ignored if algorithm is not "fit-half". Because the maximum value is virtual if lside is True and the
3652 minimum value is virtual if lside is False, the value of the maximum position (if lside=True) or
3653 minimum position (if lside=False) is not reported in the returned record.
3655 - hinges-fences: This algorithm calculates statistics by including data in a range
3656 between Q1 - f*D and Q3 + f*D, inclusive, where Q1 is the first quartile of the distribution
3657 of unmasked data, subject to any specified pixel ranges, Q3 is the third quartile, D = Q3 - Q1
3658 (the inner quartile range), and f is the user-specified fence factor. Negative values of f
3659 indicate that the full distribution is to be used (ie, the classic algorithm is used). Sufficiently
3660 large values of f will also be equivalent to using the classic algorithm. For f = 0, only data
3661 in the inner quartile range is used for computing statistics. The value of fence is silently
3662 ignored if algortihm is not "hinges-fences".
3664 - chauvenet: The idea behind this algorithm is to eliminate outliers based on a maximum z-score value.
3665 A z-score is the number of standard deviations a point is from the mean of a distribution. This
3666 method thus is meant to be used for (nearly) normal distributions. In general, this is an iterative
3667 process, with successive iterations discarding additional outliers as the remaining points become
3668 closer to forming a normal distribution. Iterating stops when no additional points lie beyond the
3669 specified zscore value, or, if zscore is negative, when Chauvenet's criterion is met (see below).
3670 The parameter maxiter can be set to a non-negative value to prematurely abort this iterative
3671 process. When verbose=True, the "N iter" column in the table that is logged represents the number
3672 of iterations that were executed.
3674 Chauvenet's criterion allows the target z-score to decrease as the number of points in the
3675 distribution decreases on subsequent iterations. Essentially, the criterion is that the probability
3676 of having one point in a normal distribution at a maximum z-score of z_max must be at least 0.5.
3677 z_max is therefore a function of (only) the number of points in the distrbution and is given by
3679 npts = 0.5/erfc(z_max/sqrt(2))
3681 where erfc() is the complementary error function. As iterating proceeds, the number of remaining
3682 points decreases as outliers are discarded, and so z_max likewise decreases. Convergence occurs when
3683 all remaining points fall within a z-score of z_max. Below is an illustrative table of z_max values
3684 and their corresponding npts values. For example, it is likely that there will be a 5-sigma "noise
3685 bump" in a perfectly noisy image with one million independent elements.
3687 z_max npts
3688 1.0 1
3689 1.5 3
3690 2.0 10
3691 2.5 40
3692 3.0 185
3693 3.5 1,074
3694 4.0 7,893
3695 4.5 73,579
3696 5.0 872,138
3697 5.5 13,165,126
3698 6.0 253,398,672
3699 6.5 6,225,098,696
3700 7.0 195,341,107,722
3702 - biweight: The biweight algorithm is a robust iterative algorithm that computes two
3703 quantities called the "location" and the "scale", which are analogous to the mean
3704 and the standard deviation. In this case, the only keys present in the returned
3705 dictionary are 'mean' (location), 'sigma' (scale), 'npts', 'min', and 'max'. The
3706 last three represent the values using the entire distribution. Note that the
3707 biweight algorithm does not support computation of quantile-like values (median,
3708 madm, q1, q3, and iqr), so setting robust=True will cause a warning message to
3709 be logged regarding that, and the computation will proceed.
3711 Important equations for the biweight algorithm are
3713 A. How to compute u_i values, which are related to the weights w_i = (1 - u_i*u_i),
3714 using the
3715 equation
3717 ``u_i = (x_i - c_bi)/(6.0*s_bi) (1)``
3719 where x_i are the data values, c_bi is the biweight location and s_bi is the
3720 biweight scale. For the initial computation of the u_i values, c_bi is set
3721 equal to the median of the distribution and s_bi is set equal to
3722 the normalized median of the absolute deviation about the median (that is the
3723 median of the absolute deviation about the median multiplied by the value of
3724 the probit function at 0.75).
3726 B. The location, c_bi, is computed from
3728 ``c_bi = sum(x_i * w_i^2)/sum(w_i^2) (2)``
3730 where only values of u_i which satisfy abs(u_i) < 1 (w_i > 0) are used in the sums.
3732 C. The scale value is computed using::
3734 . n * sum((x_i - c_bi)^2 * w_i^4)
3735 . s_bi^2 = _______________________________ (3)
3736 . p * max(1, p - 1)
3738 where n is the number of points for the entire distribution (which includes all
3739 the data, even for which abs(u_i) >= 1) and p is given by::
3741 p = abs(sum((w_i) * (5*w_i - 4)))
3743 Again, the sums include only data for which abs(u_i) < 1.
3745 The algorithm proceeds as follows.
3746 1. Compute initial u_i values (and hence w_i values) from equation (1), setting
3747 c_bi equal to the median of the distribution and s_bi equal to the normalized
3748 median of the absolute deviation about the median.
3749 2. Compute the initial value of the scale using the w_i values computed in
3750 step 1. using equation 3.
3751 3. Recompute u_i/w_i values using the most recent previous scale and location
3752 values.
3753 4. Compute the location using the u_i.w_i values from step 3 and equation (2).
3754 5. Recompute u_i/w_i values using the most recent previous scale and location
3755 values.
3756 6. Compute the new scale value using the the u_i/w_i values computed in
3757 step 5 and the value of the location computed in step 4.
3758 7. Steps 3. - 6. are repeated until convergence occurs or the maximum number of
3759 iterations (specified in the niter parameter) is reached. The convergence
3760 criterion is given by
3762 abs(s_bi - s_bi,prev) < 0.03 * sqrt(0.5/(n - 1))
3764 where s_bi,prev is the value of the scale computed in the previous iteration.
3766 In the special case where niter is specified to be negative, the faster,
3767 non-iterative algorithm proceeds as follows:
3769 1. Compute u_i/w_i values using the median for the location and the normalized
3770 median of the absolute deviation about the median as the scale
3771 2. Compute the location and scale (which can be carried out simultaneously)
3772 using the u_i/w_i values computed in step 1. The value of the location is
3773 just the median that is used in equation (3) to compute the scale
3775 NOTES ON FLUX DENSITIES AND FLUXES
3777 Fluxes and flux densities are not computed if any of the following conditions is met:
3779 1. The image does not have a direction coordinate
3780 2. The image does not have a intensity-like brightness unit. Examples of such units
3781 are Jy/beam (in which case the image must also have a beam) and K.
3782 3. There are no direction axes in the cursor axes that are used.
3783 4. If the (specified region of the) image has a non-degenerate spectral axis,
3784 and the image has a tablular spectral axis (axis with varying increments)
3785 5. Any axis that is not a direction nor a spectral axis that is included in the cursor
3786 axes is not degenerate within in the specified region
3788 Note that condition 4 may be removed in the future.
3790 In cases where none of the above conditions is met, the flux density(ies) (intensities
3791 integrated over direction planes) will be computed if any of the following conditions
3792 are met:
3794 1. The image has no spectral coordinate
3795 2. The cursor axes do not include the spectral axis
3796 3. The spectral axis in the chosen region is degenerate
3798 In the case where there is a nondegenerate spectral axis that is included in the cursor
3799 axes, the flux (flux density integrated over spectral planes) will be computed. In this
3800 case, the spectral portion of the flux unit will be the velocity unit of the spectral
3801 coordinate if it has one (eg, if the brightness unit is Jy/beam and the velocity unit is
3802 km/s, the flux will have units of Jy.km/s). If not, the spectral portion of the flux unit
3803 will be the frequency unit of the spectral axis (eg, if the brightness unit is K and the
3804 frequency unit is Hz, the resulting flux unit will be K.arcsec2.Hz).
3806 In both cases of flux density or flux being computed, the resulting numerical value is
3807 assigned to the "flux" key in the output dictionary.
3809 If the image has units of Jy/beam, the flux density is just the mean intensity multiplied
3810 by the number of beam areas included in the region. The beam area is defined as the volume
3811 of the elliptical Gaussian defined by the synthesized beam, divided by the maximum of
3812 that function, which is equivalent to
3814 pi/(4*ln(2)) * major * minor
3816 where ln() is the natural logarithm and major and minor are the major and minor FWHM axes
3817 of the beam, respectively.
3818 """
3819 return self._swigobj.statistics(axes, region, mask, includepix, excludepix, list, force, disk, robust, verbose, stretch, logfile, append, algorithm, fence, center, lside, zscore, maxiter, clmethod, niter)
3821 def twopointcorrelation(self, outfile='', region={ }, mask='', axes=[ int(-1) ], method='structurefunction', overwrite=False, stretch=False):
3822 """This function computes
3823 two-point auto-correlation functions from an image.
3825 By default, the auto-correlation function is computed for the Sky axes.
3826 If there is no sky in the image, then the first two axes are used.
3827 Otherwise you can specify which axes the auto-correlation function lags
3828 are computed over with the {stfaf axes} argument (must be of length 2).
3830 Presently, only the Structure Function is implemented. This is defined as :
3832 begin{displaymath}
3833 S(lx,ly) = < (data(i,j) - data(i+lx,j+ly))^2 >
3834 end{displaymath}
3836 where $lx, ly$ are integer lags in the x (0-axis) and y (1-axis)
3837 directions. The ensemble average is over all the values at the same
3838 lag pair. This process is extremely compute intensive and so you may
3839 have to be patient.
3841 In an auto-correlation function image there are some symmetries. The
3842 first and third quadrants are symmetric, and the second and fourth are
3843 symmetric. So in principle, all the information is in the top or bottom
3844 half of the image. We just write it all out to look nice. The long
3845 lags don't have a lot of contributing values of course.
3846 """
3847 return self._swigobj.twopointcorrelation(outfile, region, mask, axes, method, overwrite, stretch)
3849 def subimage(self, outfile='', region='', mask='', dropdeg=False, overwrite=False, list=True, stretch=False, wantreturn=True, keepaxes=[ ]):
3850 """This function copies all or part of the image to another on-the-fly Image tool.
3851 Both float and complex valued images are supported.
3853 If {stfaf outfile} is given, the subimage is written to the specified
3854 disk file. If {stfaf outfile} is unset, the returned Image tool actually
3855 references the input image file (i.e. that associated with the Image
3856 tool to which you are applying this function). So if you deleted the
3857 input image disk file, it would render this tool useless. When you
3858 destroy this tool (with the done function)
3859 the reference connection is broken.
3861 Sometimes it is useful to drop axes of length one (degenerate axes).
3862 Use the {stfaf dropdeg} argument if you want to do this. Further control
3863 is provided via the keepaxes parameter. If dropdeg=True, you may specify
3864 a list of degenerate axes to keep in the keep axes parameter. This allows
3865 you to drop only a subset of degenerate axes. This parameter is ignored if
3866 dropdeg=False. If dropdeg=True, all degenerate axes are dropped if keepaxes
3867 is set to an empty list (this is the default behavior). Nondegenerate
3868 axes are implicitly kept, even if they are included in the keepaxes list.
3870 The output mask is the combination (logical OR) of the default input
3871 pixelmask (if any) and the OTF mask. Any other input pixelmasks
3872 will not be copied. Use function maskhandler if you
3873 need to copy other masks too.
3875 If the mask has fewer dimensions than the image and if the shape
3876 of the dimensions the mask and image have in common are the same,
3877 the mask will automatically have the missing dimensions added so
3878 it conforms to the image.
3880 If stretch is true and if the number of mask dimensions is less than
3881 or equal to the number of image dimensions and some axes in the
3882 mask are degenerate while the corresponding axes in the image are not,
3883 the mask will be stetched in the degenerate dimensions. For example,
3884 if the input image has shape [100, 200, 10] and the input
3885 mask has shape [100, 200, 1] and stretch is true, the mask will be
3886 stretched along the third dimension to shape [100, 200, 10]. However if
3887 the mask is shape [100, 200, 2], stretching is not possible and an
3888 error will result.
3889 """
3890 return _wrap_image(swig_object=self._swigobj.subimage(outfile, region, mask, dropdeg, overwrite, list, stretch, wantreturn, keepaxes))
3892 def summary(self, doppler='RADIO', list=True, pixelorder=True, verbose=False):
3893 """This function summarizes various metadata such as shape, Coordinate System,
3894 restoring beams, and masks.
3896 If called without any arguments, this function displays a summary of the image
3897 metadata to the logger; where appropriate, values will be formatted nicely (e.g.
3898 HH:MM:SS.SS for the reference value of RA axes).
3900 For spectral axes, the information is listed in both the velocity and frequency
3901 domains. The doppler parameter allows one to specify what velocity doppler
3902 convention it is listed in. Supported values are: "radio", "optical", and
3903 "true". Alternative values are "z" for "optical", and "beta" or "relativistic"
3904 for true. The default is "radio". The definitions are
3906 begin{itemize}
3907 item radio: $1 - F$
3908 item optical: $-1 + 1/F$
3909 item true: $(1 - F^2)/(1 + F^2)$
3910 end{itemize}
3911 where $F = nu/nu_0$ and $nu_0$ is the rest frequency. If the rest frequency
3912 has not been set in your image, you can set it via a coordinate system (cs) tool
3913 using the setrestfrequency() method().
3915 The keys of the returned dictionary are
3917 begin{itemize}
3918 item ndim: Dimension of the image.
3919 item shape: Length of each axis in the image.
3920 item tileshape: Shape of the chunk which is most efficient for I/O.
3921 item axisnames: Name of each axis.
3922 item refpix: Reference pixel for each axis (0-relative)
3923 item refval: Reference value for each axis.
3924 item incr: Increment for each axis.
3925 item axisunits: Unit name for each axis.
3926 item unit: Brightness units for the pixels.
3927 item hasmask: True if the image has a mask.
3928 item defaultmask: The name of the mask which is applied by default.
3929 item masks: The names of all the masks stored in this image.
3930 item restoringbeam: The restoring beam(s) if present.
3931 item imagetype: The image type.
3932 end{itemize}
3934 For an image with multiple beams, the restoringbeam field is a dictionary of
3935 dictionaries with keys of names "*" followed by the channel number, if the image
3936 has a spectral coordinate, or the polarization number if it does not. That is,
3937 the keys have names "*0", "*1", "*2", etc. If the image has both a spectral and
3938 a polarization coordinate, each of these dictionaries is a dictionary with keys
3939 of the same form which range from 0 to the number of polarizations minus 1;
3940 "*0", "*1", ... The dictionaries pointed to by the channel and/or polarization
3941 number contain information for the beam at that position.
3943 If the list parameter is set to False, then the summary will not be written to
3944 the logger. The return value of the method, in the "header" field is a vector
3945 string containing the formatted output that would have been logged in the
3946 list=True case.
3948 If verbose is True and the image contains multiple beams, the formatted output,
3949 whether it is written to the logger or placed in the output record, will have
3950 information on every beam in the dataset. If verbose=False and the image has
3951 multiple beams, only a summary of beams for each polarization is listed. In this
3952 case, the beams with the maximum area, the minimum area, and the median area for
3953 each polarization are listed. However, all the beams can still be found in the
3954 restoringbeam field of the returned dictionary. If the image does not have
3955 multiple beams, verbose is not used.
3956 """
3957 return self._swigobj.summary(doppler, list, pixelorder, verbose)
3959 def tofits(self, outfile='', velocity=False, optical=True, bitpix=int(-32), minpix=float(1), maxpix=float(-1), region={ }, mask='', overwrite=False, dropdeg=False, deglast=False, dropstokes=False, stokeslast=True, wavelength=False, airwavelength=False, stretch=False, history=True):
3960 """This function converts the image into a fits file.
3963 If the image has a rest frequency associated with it, it will always
3964 write velocity information into the fits file. By default the
3965 frequency information will be primary as it is the internal native format.
3966 If you select {stfaf velocity=T} then by default
3967 the velocity is written in the optical convention, but if {stfaf
3968 optical=F} it will use the radio convention instead.
3969 Alternatively, if you use {stfaf velocity=F} and {stfaf wavelength=T},
3970 the spectral axis will be written in wavelength.
3972 The fits definition demands equal increment pixels. Therefore, if you
3973 write wavelength or optical velocity information as primary, the increment
3974 is computed at the spectral reference pixel.
3975 If the bandwidth is large, this may incur non-negligible coordinate
3976 calculation errors far from the reference pixel if the spectral
3977 bins are not originally equidistant in wavelength.
3978 Images generated by the CASA clean task have spectral axes which
3979 are always equidistant in frequency.
3981 By default the image is written as a floating point fits file
3982 ({stfaf bitpix= -32}). Under rare circumstances you might want to
3983 save space and write it as scaled 16 bit integers ({stfaf bitpix =
3984 16}). You can have {stff tofits} calculate the scaling factors by
3985 using the default {stfaf minpix} and {stfaf maxpix}. If you set
3986 {stfaf minpix} and {stfaf maxpix}, values outside of that range will
3987 be truncated. This can be useful if all of the fits images dynamic
3988 range is being used by a few high or low values and you are not
3989 interested in preserving those values exactly. Besides the factor of
3990 two space savings you get by using 16 instead of 32 bits, integer
3991 images usually also compress well (for example, with the standard GNU
3992 software facility {tt gzip}).
3994 If the specified region extends beyond the image, it is truncated.
3996 The output mask is the combination (logical OR) of the default input
3997 pixelmask (if any) and the OTF mask.
3999 Sometimes it is useful to drop axes of length one (degenerate axes)
4000 because not all FITS readers can handle them. Use the {stfaf dropdeg}
4001 argument if you want to do this.
4002 If you want to specifically only drop a degenerate Stokes axis, use the {stfaf dropstokes}
4003 argument.
4005 If you want to place degenerate axes last in the FITS header,
4006 use the {stfaf deglast} argument.
4007 If you want to make sure that the Stokes axis is placed last in the FITS header,
4008 use the {stfaf stokeslast} argument.
4009 """
4010 return self._swigobj.tofits(outfile, velocity, optical, bitpix, minpix, maxpix, region, mask, overwrite, dropdeg, deglast, dropstokes, stokeslast, wavelength, airwavelength, stretch, history)
4012 def torecord(self):
4013 """You can convert an associated image to a record for manipulation or passing it
4014 to inputs of other methods of other tools. This method and fromrecord() are used
4015 for serialization and deserialization.
4016 """
4017 return self._swigobj.torecord()
4019 def type(self):
4020 """This function returns the string 'image'. It can be used in
4021 a script to make sure this variable is an Image
4022 tool.
4023 """
4024 return self._swigobj.type()
4026 def topixel(self, value=[ ]):
4027 """This method converts from world to pixel coordinates. The world coordinate can
4028 be provided in many formats (numeric, string, quantum etc.) via the value
4029 parameter. These match the output of the toworld() method.
4031 This function is just a wrapper for the coordsys tool method of the same name,
4032 so see that documentation for a description and more examples.
4033 """
4034 return self._swigobj.topixel(value)
4036 def toworld(self, value=[ ], format='n', dovelocity=True):
4037 """This method converts between pixel and world coordinates. A variety of return
4038 formats is supported. If format='n', numerical values are returned.
4039 If format='q', values formatted as quantities are returned. If format='s',
4040 values formatted as strings are returned. If format='m', values formatted as
4041 measures are returned. If format='m', one can choose to have the corresponding
4042 velcocites of an extant spectral coordinate computed as well by specifyting
4043 dovelocity=True (dovelocity is ignored if format is not equal to 'm' or if the
4044 image does not have a spectral coordinate).
4045 """
4046 return self._swigobj.toworld(value, format, dovelocity)
4048 def unlock(self):
4049 """This function releases any lock set on the imagefile (and also flushes
4050 any outstanding I/O to disk). It is not of general user interest. It
4051 can be useful in scripts when a file is being shared between more than
4052 one process. See also functions lock and
4053 haslock.
4054 """
4055 return self._swigobj.unlock()
4057 def newimagefromarray(self, outfile='', pixels=[ ], csys={ }, linear=False, overwrite=False, log=True, type='f'):
4058 """This application converts a numpy array of any size into a CASA image.
4060 If outfile is specified, the image is written to the specified
4061 (persistent) disk file. If outfile is unset, the returned image tool
4062 is associated with a temporary image. This temporary image may be in
4063 memory or on disk, depending on its size. In this case, when the
4064 close() or done() method is called on the returned image tool, the
4065 associated temporary image is deleted.
4067 The type parameter controls the data type/precision of the pixel values of the
4068 created image. 'f' indicates that float precision point (32 bit precision) pixel
4069 values should be writted. 'd' indicates that double precision (64 bit precision)
4070 pixel values should be written. If the input array has complex (as opposed to
4071 real) values, then complex pixel values, with each of the real and imaginary
4072 parts having the specified precision, will be written. Array values will be cast
4073 automatically to the specified precision, so that the precision of the input
4074 array values may be increased, decreased, or unchanged depending on the input
4075 array type.
4077 The coordinate system, provided as a a dictionary (use eg, cs.torecord() to do
4078 that), is optional. If specified, it must have the same number of dimensions
4079 as the pixels array. Call the naxes() method on the coordinate system tool to
4080 see how many dimensions the coordinate system has. A coordinate system can be
4081 created from scratch using the coordinate system (cs) tool and methods therein,
4082 but often users prefer to use a coordinate system from an already existing image.
4083 This can be gotten using ia.coordsys() which returns a coordinate system tool.
4084 A torecord() call on that tool will result in a python dictionary describing
4085 the coordinate system which is the necessary format for the csys input parameter
4086 of ia.fromarray().
4088 If csys is not specified, a default coordinate system will be created. If
4089 linear=False (the default), the created coordinate system will have standard
4090 RA/DEC/Stokes/Spectral Coordinate axes depending upon the shape of the pixels
4091 array (Stokes axis must be no longer than 4 pixels and the spectral axis may
4092 precede the Stokes axis if eg, shape=[64,64,32,4]. Extra dimensions are given
4093 linear coordinates. If linear=True, then all the resulting coordinates
4094 are linear with the axes represent lengths. In this case each axis will have a
4095 value of 0.0 at its center pixel. The increment of each axis will be 1.0 km.
4096 """
4097 return _wrap_image(swig_object=self._swigobj.newimagefromarray(outfile, pixels, csys, linear, overwrite, log, type))
4099 def newimagefromfits(self, outfile='', infile='', whichrep=int(0), whichhdu=int(0), zeroblanks=False, overwrite=False):
4100 """This function is used to convert a FITS disk image file (Float,
4101 Double, Short, Long are supported) to an casa imagefile. If
4102 {stfaf outfile} is given, the image is written to the specified disk
4103 file. If {stfaf outfile} is unset, the on-the-fly Image tool
4104 returned by this function is associated with a temporary image. This
4105 temporary image may be in memory or on disk, depending on its size.
4106 When you destroy the on-the-fly Image tool (with the done function) this
4107 temporary image is deleted.
4109 This function reads from the FITS primary array (when the image is at
4110 the beginning of the FITS file; {stfaf whichhdu=0}), or an image
4111 extension (when the image is elsewhere in the FITS file, {stfaf
4112 whichhdu $>$ 0}).
4114 By default, any blanked pixels will be converted to a mask value which
4115 is false, and a pixel value that is NaN. If you set {stfaf
4116 zeroblanks=T} then the pixel value will be zero rather than NaN. The
4117 mask will still be set to false. See the function
4118 replacemaskedpixels if you
4119 need to replace masked pixel values after you have created the image.
4120 """
4121 return _wrap_image(swig_object=self._swigobj.newimagefromfits(outfile, infile, whichrep, whichhdu, zeroblanks, overwrite))
4123 def newimagefromimage(self, infile='', outfile='', region={ }, mask='', dropdeg=False, overwrite=False):
4124 """This function applies a region to a disk imagefile, creates a new
4125 imagefile containing the (sub)image, and associates a new imagetool
4126 with it.
4128 The input disk image file may be in native casa, fits (Float,
4129 Double, Short, Long are supported), or Miriad format. Look
4130 htmlref{here}{IMAGES:FOREIGNIMAGES} for more information on foreign
4131 images.
4133 If {stfaf outfile} is given, the (sub)image is written to the specified
4134 disk file.
4136 If {stfaf outfile} is unset, the Image tool actually references the
4137 input image file. So if you deleted the input image disk file, it
4138 would render this tool useless. When you destroy this on-the-fly
4139 tool (with the done
4140 function) the reference connection is broken.
4142 Sometimes it is useful to drop axes of length one (degenerate axes).
4143 Use the {stfaf dropdeg} argument if you want to do this.
4145 The output mask is the combination (logical OR) of the default input
4146 pixelmask (if any) and the OTF mask. Any other input pixelmasks
4147 will not be copied. Use function
4148 maskhandler if you need to copy other
4149 masks too.
4151 See also the subimage function.
4152 """
4153 return _wrap_image(swig_object=self._swigobj.newimagefromimage(infile, outfile, region, mask, dropdeg, overwrite))
4155 def newimagefromshape(self, outfile='', shape=[ int(0) ], csys={ }, linear=False, overwrite=False, log=True, type='f'):
4156 """This function creates a CASA image with the specified shape. It is similar to
4157 ia.fromshape(), but instead returns a new image analysis tool attached to the
4158 new image, rather than attaching the new image to the current tool. All the
4159 pixel values in the image are set to 0. One may create an image with float
4160 precision pixels (type='f'), complex float precision pixels (type='c'), double
4161 precision pixels (type='d'), or complex double precision pixels ('cd'). To use a
4162 numpy array of values to create an image, use ia.(newimage)fromarray(). To make
4163 a 2-D image from a packaged FITS file, use ia.maketestimage().
4165 If outfile is given, the image is written to the specified disk file. If
4166 outfile is unset, the image analysis tool is associated with a temporary image.
4167 This temporary image may be in memory or on disk, depending on its size. When
4168 you close the image analysis tool (with the ia.close() method, the temporary
4169 image is deleted.
4171 The coordinate system, provided as a coordinate system tool record, is optional.
4172 If provided, it must be dimensionally consistent with the specified shape.
4174 If the coordinate system is not provided, a default coordinate system will be
4175 created. If linear=False (the default), then it is a
4176 standard RA/DEC/Stokes/Spectral coordinate system depending exactly upon the
4177 shape (the Stokes axis must be no longer than 4 pixels and spectral axis may
4178 occur prior to the Stokes axis if eg, shape=[64,64,32,4]. Extra dimensions are
4179 given linear coordinates. If linear=True, then the coordinate system will have
4180 linear coordinates.
4181 """
4182 return _wrap_image(swig_object=self._swigobj.newimagefromshape(outfile, shape, csys, linear, overwrite, log, type))
4184 def pbcor(self, pbimage='', outfile='', overwrite=False, box='', region={ }, chans='', stokes='', mask='', mode='divide', cutoff=float(-1.0), stretch=False):
4185 """Correct an image for primary beam attenuation using an image of the primary beam pattern.
4186 The primary beam pattern can be provided as an image, in which case 1. it must have the same
4187 shape as the input image and its coordinate system must be the same, or 2. it must
4188 be a 2-D image in which case its coordinate system must consist of a (2-D) direction
4189 coordinate which is the same as the direction coordinate in the input image and
4190 its direction plane must be the same shape as that of the input image. Alternatively,
4191 pbimage can be an array of pixel values in which case the same dimensionality and
4192 shape constraints apply.
4193 An image tool referencing the corrected image is returned. The corrected image will also
4194 be written to disk if outfile is not empty (and overwrite=True if outfile already exists).
4195 One can choose between dividing the image by the primary beam pattern (mode="divide") or
4196 multiplying the image by the primary beam pattern (mode="multiply"). One can also choose
4197 to specify a cutoff limit for the primary beam pattern. For mode="divide", for all pixels
4198 below this cutoff in the primary beam pattern, the output image will be masked. In the
4199 case of mode="multiply", all pixels in the output will be masked corresponding to pixels
4200 with values greater than the cutoff in the primary beam pattern. A negative value for
4201 cutoff means that no cutoff will be applied, which is the default.
4203 """
4204 return _wrap_image(swig_object=self._swigobj.pbcor(pbimage, outfile, overwrite, box, region, chans, stokes, mask, mode, cutoff, stretch))
4206 def pixeltype(self):
4207 """This application returns the data type of the pixels of the attached image as a string.
4208 The possible values are: "float" which indicates real valued, floating point, 32 bit pixel
4209 values, "complex" which indicates complex valued, floating point, 32 bit (for each of the
4210 real and imaginary parts) pixel values, "double" which indicates real valued, floating
4211 point, 64 bit pixel values, and "dcomplex" which indicates complex valued, floating point,
4212 64 bit (for each of the real and imaginary parts) pixel values.
4214 """
4215 return self._swigobj.pixeltype()
4217 def pv(self, outfile='', start=[ ], end=[ ], center=[ ], length=[ ], pa=[ ], width=int(1), unit='arcsec', overwrite=False, region={ }, chans='', stokes='', mask='', stretch=False, wantreturn=True):
4218 """Create a position-velocity image by specifying either two points between which a slice is taken in the direction
4219 coordinate or a center, position angle, and length describing the slice. The spectral extent of the resulting image
4220 will be that provided by the region specification or the entire spectral range of the input image if no region is
4221 specified. One may not specify a region in direction space; that is accomplished by specifying the start and end
4222 points or the center, length, and position angle of the slice. The parameters start and end may be specified as two
4223 element arrays of numerical values, in which case these values will be interpreted as pixel locations in the input
4224 image. Alternatively, they may be expressed as arrays of two strings each representing the direction. These strings
4225 can either represent quantities (eg ["40.5deg", "0.5rad") or be sexigesimal format (eg ["14h20m20.5s","-30d45m25.4s"],
4226 ["14:20:20.5s","-30.45.25.4"]). In addition, they may be expressed as a single string containing the longitude and
4227 latitude values and optionally a reference frame value, eg "J2000 14:20:20.5s -30.45.25.4". The center parameter is
4228 specified in the same way. The length parameter may be specified as a single numerical value, in which case it is
4229 interpreted as the length in pixels, or a valid quantity, in which case it must have units conformant with the direction
4230 axes units. The pa (position angle) parameter must be specified as a valid quantity with angular units. The position
4231 angle is interpreted in the usual astronomical sense; ie measured from north through east. Either start/end or
4232 center/pa/length must be specified; if a parameter from one of these sets is specified, a parameter from the other set may
4233 not be specified. In either case, the end points of the segment must fail within the input image, and they both must be at
4234 least 2 pixels from the edge of the input image to facilite rotation (see below).
4236 One may specify a width, which is the number of pixels centered along and perpendicular
4237 to the direction slice that are used for averaging along the slice. The width may be specified as an integer, in which
4238 case it must be positive and odd. Alternatively, it may be specified as a valid quantity string (eg, "4arcsec") or
4239 quantity record (eg qa.quantity("4arcsec"). In this case, units must be conformant to the direction axes units (usually
4240 angular units) and the specified quantity will be rounded up, if necessary, to the next highest equivalent odd integer number
4241 of pixels. The default value of 1 represents no averaging.
4242 A value of 3 means average one pixel on each side of the slice and the pixel on the slice.
4243 Note that this width is applied to pixels in the image after it has been rotated (see below for a description
4244 of the algorithm used). The end points of the specified segment must fail within the input
4245 image, and they both must be at least 2 pixels from the edge of the input image to facilite rotation (see below).
4247 One may specify the unit for the angular offset axis.
4249 A true value for the wantreturn parameter indicates that an image analysis tool attached to the created
4250 image should be returned. Nothing is returned if wantreturn is false, but then outfile should be specified
4251 (unless perhaps you are debugging).
4253 Internally, the image is first rotated, padding first if necessary to include relevant pixels that would otherwise
4254 be excluded by the rotation operation, so that the slice is horizontal, with the starting pixel left of the
4255 ending pixel. Then, the pixels within the specified width of the slice are averaged and the resulting image is
4256 written and/or returned. The output image has a linear coordinate in place of the direction coordinate of the
4257 input image, and the corresponding axis represents angular offset with the center pixel having a value of 0.
4259 The equivalent coordinate system, with a (usually) rotated direction coordinate (eg, RA and Dec) is written
4260 to the output image as a table record. It can be retrieved using the table tool as shown in the example below.
4262 """
4263 return _wrap_image(swig_object=self._swigobj.pv(outfile, start, end, center, length, pa, width, unit, overwrite, region, chans, stokes, mask, stretch, wantreturn))
4265 def makearray(self, v=float(0.0), shape=[ int(0) ]):
4266 """This function takes two arguments. The first argument is the initial
4267 value for the new array. The second is a vector giving the lengths of
4268 the dimensions of the array.
4269 """
4270 return self._swigobj.makearray(v, shape)
4272 def isconform(self, other):
4273 """Returns True if the shape, coordinate system, and axes order of the specified image
4274 matches the current image.
4275 """
4276 return self._swigobj.isconform(other)