Coverage for /wheeldirectory/casa-6.7.0-12-py3.10.el8/lib/py/lib/python3.10/site-packages/casatasks/private/task_importfitsidi.py: 89%
108 statements
« prev ^ index » next coverage.py v7.6.4, created at 2024-11-01 07:19 +0000
1import os
2import shutil
3import numpy as np
5from casatools import ms, table
6from casatasks import casalog
7from .mstools import write_history
9def importfitsidi(fitsidifile,vis,constobsid=None,scanreindexgap_s=None,specframe=None):
10 """Convert FITS-IDI visibility file into a CASA visibility file (MS).
12 Keyword arguments:
13 fitsidifile -- Name(s) of input FITS IDI file(s)
14 default: None; example='3C273XC1.IDI' or ['3C273XC1.IDI1', '3C273XC1.IDI2']
15 vis -- Name of output visibility file (MS)
16 default: None; example: vis='3C273XC1.ms'
18 constobsid -- If True a constant obs id == 0 of is given to all input files
19 default = False (new obs id for each input file)
21 scanreindexgap_s -- if > 0., a new scan is started whenever the gap between two
22 integrations is > the given value (seconds) or when a new field starts
23 or when the ARRAY_ID changes.
24 default = 0. (no reindexing)
26 specframe -- this frame will be used to set the spectral reference frame
27 for all spectral windows in the output MS
28 default = GEO (geocentric), other options: TOPO, LSRK, BARY
29 NOTE: if specframe is set to TOPO, the reference location will be taken from
30 the Observatories table in the CASA data repository for the given name of
31 the observatory. You can edit that table and add new rows.
33 """
35 #Python script
37 try:
38 casalog.origin('importfitsidi')
39 casalog.post("")
40 myms = ms()
41 mytb = table()
43 if type(specframe)==str and not specframe=='':
44 myspecframe=specframe.upper()
45 else:
46 myspecframe='GEO'
48 refframes = {'REST': 0, 'LSRK': 1, 'LSRD': 2, 'BARY': 3, 'GEO': 4, 'TOPO': 5}
49 if myspecframe not in refframes:
50 raise ValueError('Value '+myspecframe+' of parameter specframe invalid. Possible values are REST, LSRK, LSRD, BARY, GEO, TOPO')
52 if(type(fitsidifile)==str):
53 casalog.post('### Reading file '+fitsidifile, 'INFO')
54 myms.fromfitsidi(vis,fitsidifile)
55 myms.close()
56 elif(type(fitsidifile)==list):
57 clist = fitsidifile
58 casalog.post('### Reading file '+clist[0], 'INFO')
59 myms.fromfitsidi(vis,clist[0])
60 myms.close()
61 clist.pop(0)
62 tname = vis+'_importfitsidi_tmp_'
63 shutil.rmtree(tname, ignore_errors=True)
64 for fidifile in clist:
65 casalog.post('### Reading file '+fidifile, 'INFO')
66 myms.fromfitsidi(tname,fidifile)
67 myms.close()
68 myms.open(vis, nomodify=False)
69 myms.concatenate(msfile=tname, freqtol='', dirtol='')
70 myms.close()
71 shutil.rmtree(tname, ignore_errors=True)
72 else:
73 raise ValueError('Parameter fitsidifile should be of type str or list')
75 if (constobsid):
76 mytb.open(vis+'/OBSERVATION', nomodify=False)
77 nobs = mytb.nrows()
78 cando = True
79 if nobs>1:
80 casalog.post('Trying to keep obsid constant == 0 for all input files', 'INFO')
81 # check if all observations are from the same telescope; if not warn and leave as is
82 tels = mytb.getcol('TELESCOPE_NAME')
83 for i in range(1,nobs):
84 if tels[i]!=tels[0]:
85 cando = False
87 if cando:
88 # get min and max time and write them into the first row;
89 casalog.post('Adjusting OBSERVATION table', 'INFO')
90 timeranges = mytb.getcol('TIME_RANGE')
91 newmin = min(timeranges[0])
92 newmax = max(timeranges[1])
93 mytb.putcell('TIME_RANGE', 0, [newmin,newmax])
94 # delete the other rows
95 mytb.removerows(list(range(1,nobs)))
96 else:
97 casalog.post('The input files stem from different telescopes. Need to give different obs id.', 'WARN')
98 mytb.close()
100 if cando:
101 # give the same obs id == 0 to the entire output MS
102 casalog.post('Setting observation ID of all integrations to 0', 'INFO')
103 mytb.open(vis, nomodify=False)
104 for i in range(0, mytb.nrows()):
105 mytb.putcell('OBSERVATION_ID', i, 0)
106 mytb.close()
109 else: # don't want constant obs id
110 if(type(fitsidifile)==list and len(fitsidifile)>1):
111 casalog.post('Incrementing observation ID for each input file ...', 'INFO')
113 if (scanreindexgap_s > 0.):
114 # reindex the scan column
115 mytb.open(vis, nomodify=False)
116 times = mytb.getcol('TIME')
117 fields = mytb.getcol('FIELD_ID')
118 arrayids = mytb.getcol('ARRAY_ID')
119 scannumbers = mytb.getcol('SCAN_NUMBER')
121 timesorted = np.argsort(np.array(times))
123 scannumber = 0
124 scannumber_field = len(fields) * [0]
125 prevtime = len(fields) * [0]
126 prevarrayid = arrayids[timesorted[0]]
128 for i in range(0,mytb.nrows()):
129 ii = timesorted[i]
130 timenow = times[ii]
131 fieldnow = fields[ii]
132 arrayidnow = arrayids[ii]
133 if (timenow-prevtime[fieldnow] > scanreindexgap_s) \
134 or (arrayidnow != prevarrayid):
135 scannumber += 1
136 scannumber_field[fieldnow] = scannumber
137 casalog.post("Starting new scan "+str(scannumber)+" at "+str(timenow)\
138 +", field "+str(fieldnow)+", array_id "+str(arrayidnow), 'INFO')
139 scannumbers[ii] = scannumber_field[fieldnow]
140 prevtime[fieldnow] = timenow
141 prevarrayid = arrayidnow
143 mytb.putcol('SCAN_NUMBER', scannumbers)
144 mytb.close()
146 if myspecframe in refframes:
147 casalog.post('Setting reference frame for all spectral windows to '+myspecframe, 'INFO')
148 if myspecframe == 'TOPO':
149 casalog.post('NOTE: reference position for TOPO frame will be the observatory location', 'WARN')
150 mytb.open(vis+'/SPECTRAL_WINDOW', nomodify=False)
151 refcol = mytb.getcol('MEAS_FREQ_REF')
152 refcol = [refframes[myspecframe]]*len(refcol)
153 mytb.putcol('MEAS_FREQ_REF', refcol)
154 mytb.close()
156 # write history
157 try:
158 param_names = importfitsidi.__code__.co_varnames[:importfitsidi.__code__.co_argcount]
159 vars = locals( )
160 param_vals = [vars[p] for p in param_names]
161 write_history(myms, vis, 'importfitsidi', param_names, param_vals, casalog)
163 except Exception as instance:
164 casalog.post("*** Error \'%s\' updating HISTORY" % instance, 'WARN')
166 finally:
167 shutil.rmtree(vis+'_importfitsidi_tmp_', ignore_errors=True)