casa
$Rev:20696$
|
00001 ########################################################################3 00002 # task_immoments.py 00003 # 00004 # 00005 # Copyright (C) 2008, 2009 00006 # Associated Universities, Inc. Washington DC, USA. 00007 # 00008 # This script is free software; you can redistribute it and/or modify it 00009 # under the terms of the GNU Library General Public License as published by 00010 # the Free Software Foundation; either version 2 of the License, or (at your 00011 # option) any later version. 00012 # 00013 # This library is distributed in the hope that it will be useful, but WITHOUT 00014 # ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or 00015 # FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public 00016 # License for more details. 00017 # 00018 # You should have received a copy of the GNU Library General Public License 00019 # along with this library; if not, write to the Free Software Foundation, 00020 # Inc., 675 Massachusetts Ave, Cambridge, MA 02139, USA. 00021 # 00022 # Correspondence concerning AIPS++ should be adressed as follows: 00023 # Internet email: aips2-request@nrao.edu. 00024 # Postal address: AIPS++ Project Office 00025 # National Radio Astronomy Observatory 00026 # 520 Edgemont Road 00027 # Charlottesville, VA 22903-2475 USA 00028 # 00029 # <author> 00030 # Shannon Jaeger (University of Calgary) 00031 # </author> 00032 # 00033 # <summary> 00034 # CASA task for finding moments along a specified axis of a 00035 # multi-dimentional CASA image. 00036 # contents 00037 # </summary> 00038 # 00039 # <reviewed reviwer="" date="" tests="" demos=""> 00040 # </reviewed 00041 00042 # <etymology> 00043 # immoments stands for image momemnts 00044 # </etymology> 00045 # 00046 # <synopsis> 00047 # task_immoments.py is a Python script providing an easy to use task 00048 # for generating momements along a user specified axis on a CASA image. 00049 # This is a time-honoured spectral-line analysis technique for discovering 00050 # spectral line information. 00051 # 00052 # In this task, moment, refers to collapsing an axis of the image, 00053 # the moment axis, to a single pixel. 00054 # 00055 # The various moments that can be calculated are described in detail 00056 # at http://casa.nrao.edu/docs/casaref/image.moments.html#x59-590001.1.1 00057 # 00058 # </synopsis> 00059 # 00060 # <example> 00061 # <srcblock> 00062 # # The following code snippet find the 1-moments, intensity-weighted 00063 # # coordinate, often used for finding velocity fields. 00064 # immoments( 'axis='spec', imagename='myimage', moment=1, outfile='velocityfields' ) 00065 # 00066 # # Finding the spectral mean, -1 moment, on a specified portion of the image. 00067 # # The region used is defined by the box and stokes parameters 00068 # immoments( imagename='myimage', axis='spec', stokes='I', box=[55,12,97,32], moment=-1 ) 00069 # 00070 # # The following example uses a second file to use as a mask. 00071 # # The 0-moments, integrated values, are created on clean.im, but 00072 # # the mask is based on all the data in calibrated.im, all values 00073 # # about the threshold 0.5 are used to create the mask. 00074 # 00075 # immoments( 'clean.image', axis='spec', mask='calibrated.im>0.5', outfile='mom_withmask.im' ) 00076 # </srblock> 00077 # </example> 00078 # 00079 # <motivation> 00080 # To provide a user-friendly method of calculating image moments. 00081 # </motivation> 00082 # 00083 # <todo> 00084 # </todo> 00085 00086 import os 00087 from taskinit import * 00088 00089 def immoments( 00090 imagename, moments, axis, region, box, chans, stokes, 00091 mask, includepix, excludepix, outfile, stretch 00092 ): 00093 00094 retValue=None 00095 casalog.origin('immoments') 00096 00097 # First check to see if the output file exists. If it 00098 # does then we abort. CASA doesn't allow files to be 00099 # over-written, just a policy. 00100 if ( len( outfile ) > 0 and os.path.exists( outfile ) ): 00101 raise Exception, 'Output file, '+outfile+\ 00102 ' exists. immoment can not proceed, please\n'\ 00103 'remove it or change the output file name.' 00104 elif ( len( outfile ) < 1 ): 00105 casalog.post( "The outfile paramter is empty, consequently the" \ 00106 +" moment image will NOT be\nsaved on disk," \ 00107 +" but an image tool (ia) will be returned and if the" \ 00108 +" returned value\nis saved then you can used in" \ 00109 +" the same way the image tool (ia). can", 'WARN' ) 00110 00111 _myia = iatool() 00112 try: 00113 # Translate the string value into an index value. 00114 axis=_immoments_parse_axis( imagename, axis ) 00115 casalog.post( 'Axis information for '+imagename+' is: '+str(axis),\ 00116 'DEBUG2' ) 00117 _myia.open(imagename) 00118 00119 reg = rg.frombcs(csys=_myia.coordsys().torecord(), 00120 shape=_myia.shape(), box=box, chans=chans, stokes=stokes, 00121 stokescontrol="a", region=region 00122 ) 00123 retValue = _myia.moments( 00124 moments=moments, axis=int(axis), mask=mask, 00125 region=reg, includepix=includepix, 00126 excludepix=excludepix, outfile=outfile, drop=False, 00127 stretch=stretch 00128 ) 00129 _myia.done() 00130 return retValue 00131 except Exception, instance: 00132 _myia.done() 00133 print '*** Error ***',instance 00134 raise Exception, instance 00135 00136 00137 def _immoments_parse_axis( imagename='', axis='' ): 00138 # We already have an integer value nothing to do. 00139 if ( isinstance( axis, int ) ): 00140 return axis 00141 00142 # Find out what we have in each axis of this image. 00143 # of a particular region in the image. This function 00144 # returns the following list: 00145 # [ DirectionalIndex, Directional1Index, spectralIndex, stokesIndex ] 00146 axes=getimaxes(imagename) 00147 00148 # Default is spectral axis 00149 retValue=axes[2] 00150 00151 00152 # Ignore case and remove extra spaces 00153 axis = axis.replace( ' ', '' ) 00154 axis = axis.lower() 00155 00156 00157 if ( axis=='ra' or axis.startswith( 'long' ) ): 00158 # First Directional axis 00159 retValue=axes[0][0] 00160 elif ( axis=='dec' or axis.startswith( 'lat' ) ): 00161 # Second Directional axis 00162 retValue=axes[1][0] 00163 elif ( axis.startswith( 'spec' ) ): 00164 # Spectral axis 00165 retValue=axes[2][0] 00166 elif ( axis.startswith( 'sto') ): 00167 # Spectral axis 00168 retValue=axes[3][0] 00169 else: 00170 raise Exception, "Invalid axis specified: "+str(axis) \ 00171 + ". Expected one of ra, dec, lat, long, spec, or stokes" 00172 00173 return retValue 00174