Coverage for /wheeldirectory/casa-6.7.0-12-py3.10.el8/lib/py/lib/python3.10/site-packages/casatasks/private/task_concat.py: 7%
226 statements
« prev ^ index » next coverage.py v7.6.4, created at 2024-10-31 17:39 +0000
1import os
2import shutil
3import stat
4import time
5from math import sqrt
7from .parallel.parallel_task_helper import ParallelTaskHelper
8from .mslisthelper import check_mslist, sort_mslist
9from casatools import calibrater, quanta
10from casatools import table as tbtool
11from casatools import ms as mstool
12from casatasks import casalog
13from .mstools import write_history
15_cb = calibrater()
16_qa = quanta()
18def concat(vislist,concatvis,freqtol,dirtol,respectname,timesort,copypointing,
19 visweightscale, forcesingleephemfield):
20 """concatenate visibility datasets
21 The list of data sets given in the vis argument are chronologically concatenated
22 into an output data set in concatvis, i.e. the data sets in vis are first ordered
23 by the time of their earliest integration and then concatenated.
25 If there are fields whose direction agrees within the direction tolerance
26 (parameter dirtol), the actual direction in the resulting, merged output field
27 will be the one from the chronologically first input MS.
29 If concatvis already exists (e.g., it is the same as the first input data set),
30 then the other input data sets will be appended to the concatvis data set.
31 There is no limit to the number of input data sets.
33 If none of the input data sets have any scratch columns (model and corrected
34 columns), none are created in the concatvis. Otherwise these columns are
35 created on output and initialized to their default value (1 in model column,
36 data in corrected column) for those data with no input columns.
38 Spectral windows for each data set with the same chanelization, and within a
39 specified frequency tolerance of another data set will be combined into one
40 spectral window.
42 A field position in one data set that is within a specified direction tolerance
43 of another field position in any other data set will be combined into one
44 field. The field names need not be the same---only their position is used.
46 Each appended dataset is assigned a new observation id (provided the entries
47 in the observation table are indeed different).
49 Keyword arguments:
50 vis -- Name of input visibility files to be combined
51 default: none; example: vis = ['src2.ms','ngc5921.ms','ngc315.ms']
52 concatvis -- Name of visibility file that will contain the concatenated data
53 note: if this file exits on disk then the input files are
54 added to this file. Otherwise the new file contains
55 the concatenated data. Be careful here when concatenating to
56 an existing file.
57 default: none; example: concatvis='src2.ms'
58 example: concatvis='outvis.ms'
60 other examples:
61 concat(vis=['src2.ms','ngc5921.ms'], concatvis='src2.ms')
62 will concatenate 'ngc5921.ms' into 'src2.ms', and the original
63 src2.ms is lost
65 concat(vis=['src2.ms','ngc5921.ms'], concatvis='out.ms')
66 will concatenate 'ngc5921.ms' and 'src2.ms' into a file named
67 'out.ms'; the original 'ngc5921.ms' and 'src2.ms' are untouched.
69 concat(vis=['v1.ms','v2.ms'], concatvis = 'vall.ms')
70 then
71 concat(vis=['v3.ms','v4.ms'], concatvis = 'vall.ms')
72 vall.ms will contains v1.ms+v2.ms+v3.ms+v4.ms
74 Note: run flagmanager to save flags in the concatvis
76 freqtol -- Frequency shift tolerance for considering data to be in the same
77 spwid. The number of channels must also be the same.
78 default: '' == 1 Hz
79 example: freqtol='10MHz' will not combine spwid unless they are
80 within 10 MHz.
81 Note: This option is useful to combine spectral windows with very slight
82 frequency differences caused by Doppler tracking, for example.
84 dirtol -- Direction shift tolerance for considering data as the same field
85 default: '' == 1 mas (milliarcsec)
86 example: dirtol='1.arcsec' will not combine data for a field unless
87 their phase center differ by less than 1 arcsec. If the field names
88 are different in the input data sets, the name in the output data
89 set will be the first relevant data set in the list.
91 respectname -- If true, fields with a different name are not merged even if their
92 direction agrees (within dirtol)
93 default: False
95 timesort -- If true, the output visibility table will be sorted in time.
96 default: false. Data in order as read in.
97 example: timesort=true
98 Note: There is no constraint on data that is simultaneously observed for
99 more than one field; for example multi-source correlation of VLBA data.
101 copypointing -- Make a proper copy of the POINTING subtable (can be time consuming).
102 If False, the result is an empty POINTING table.
103 default: True
105 visweightscale -- The weights of the individual MSs will be scaled in the concatenated
106 output MS by the factors in this list. Useful for handling heterogeneous arrays.
107 Use plotms to inspect the "Wt" column as a reference for determining the scaling
108 factors. See the cookbook for more details.
109 example: [1.,3.,3.] - scale the weights of the second and third MS by a factor 3.
110 default: [] (empty list) - no scaling
112 forcesingleephemfield -- By default, concat will only merge two ephemeris fields if
113 the first ephemeris covers the time range of the second. Otherwise, two separate
114 fields with separate ephemerides are placed in the output MS.
115 In order to override this behaviour and make concat merge the non-overlapping
116 or only partially overlapping input ephemerides, the name or id of the field
117 in question needs to be placed into the list in parameter 'forcesingleephemfield'.
118 example: ['Neptune'] - will make sure that there is only one joint ephemeris for
119 field Neptune in the output MS
120 default: '' - standard treatment of all ephemeris fields
122 """
124 ###
125 #Python script
126 try:
127 casalog.origin('concat')
128 t = tbtool()
129 m = mstool()
131 #break the reference between vis and vislist as we modify vis
132 if(type(vislist)==str):
133 vis=[vislist]
134 else:
135 vis=list(vislist)
136 #dto. for concavis
137 theconcatvis = concatvis
139 # warn if there are MMSs
140 mmslist = []
141 for elvis in vis : ###Oh no Elvis does not exist Mr Bill
142 if(ParallelTaskHelper.isParallelMS(elvis)):
143 mmslist.append(elvis)
144 if len(mmslist)>0:
145 if (vis[0] == mmslist[0]):
146 casalog.post('*** The first input MS is a multi-MS to which no row can be added. Cannot proceed.', 'WARN')
147 casalog.post('*** Please use virtualconcat or convert the first input MS to a normal MS using split.', 'WARN')
148 raise RuntimeError('Cannot append to a multi-MS. Please use virtualconcat.')
150 casalog.post('*** The following input measurement sets are multi-MSs', 'INFO')
151 for mname in mmslist:
152 casalog.post('*** '+mname, 'INFO')
153 casalog.post('*** Use virtualconcat to produce a single multi-MS from several multi-MSs.', 'INFO')
157 doweightscale = False
158 if(len(visweightscale)>0):
159 if (len(visweightscale) != len(vis)):
160 raise ValueError('parameter visweightscale must have same number of elements as parameter vis')
161 for factor in visweightscale:
162 if factor<0.:
163 raise ValueError('parameter visweightscale must only contain positive numbers')
164 elif factor!=1.:
165 doweightscale=True
168 # test the consistency of the setup of the different MSs
169 casalog.post('Checking MS setup consistency ...', 'INFO')
170 try:
171 mydiff = check_mslist(vis, ignore_tables=['SORTED_TABLE', 'ASDM*'], testcontent=False)
172 except Exception as instance:
173 raise RuntimeError("*** Error \'%s\' while checking MS setup consistency" % (instance))
175 if mydiff != {}:
176 casalog.post('The setup of the input MSs is not fully consistent. The concatenation may fail', 'WARN')
177 casalog.post('and/or the affected columns may contain partially only default data.', 'WARN')
178 casalog.post(str(mydiff), 'WARN')
181 # process the input MSs in chronological order
182 casalog.post('Checking order of MS list ...', 'INFO')
183 try:
184 sortedvis, sortedtimes, sortedvisweightscale = sort_mslist(vis, visweightscale)
185 except Exception as instance:
186 raise RuntimeError("*** Error \'%s\' while sorting MSs chronologially." % (instance))
188 if((type(concatvis)!=str) or (len(concatvis.split()) < 1)):
189 raise ValueError('parameter concatvis is invalid')
191 existingconcatvis = False
192 if(vis.count(concatvis) > 0):
193 existingconcatvis = True
194 cvisindex = sortedvis.index(concatvis)
195 if not concatvis == sortedvis[0]:
196 raise RuntimeError('If concatvis is set to the name of an existing MS in vis, it must be the chronologically first.'+\
197 '\n I.e. in this case you should set concatvis to '+sortedvis[0])
198 sortedvis.pop(cvisindex)
199 if doweightscale:
200 vwscale = sortedvisweightscale[cvisindex]
201 sortedvisweightscale.pop(cvisindex)
202 sortedvisweightscale = [vwscale] + sortedvisweightscale # move the corresponding weight to the front
204 if not vis == sortedvis:
205 casalog.post('The list of input MSs is not in chronological order and needed to be sorted.' , 'INFO')
206 casalog.post('The chronological order in which the concatenation will take place is:' , 'INFO')
207 if existingconcatvis:
208 casalog.post(' MJD '+str(_qa.splitdate(_qa.quantity(sortedtimes[0],'s'))['mjd'])+': '+concatvis, 'INFO')
209 for name in sortedvis:
210 casalog.post(' MJD '+str(_qa.splitdate(_qa.quantity(sortedtimes[sortedvis.index(name)],'s'))['mjd'])+': '+name, 'INFO')
211 if doweightscale:
212 casalog.post('In this new order, the weights are:'+str(sortedvisweightscale) , 'INFO')
214 # replace the original vis and visweightscale by the sorted ones (with concatvis removed if it exists)
215 vis = sortedvis
216 visweightscale = sortedvisweightscale
218 if(os.path.exists(concatvis)):
219 casalog.post('Will be concatenating into the existing ms '+concatvis , 'WARN')
220 if doweightscale and not existingconcatvis:
221 visweightscale = [1.]+visweightscale # set the weight for this existing MS to 1.
222 casalog.post('The weights for this existing MS will be left unchanged.' , 'WARN')
223 else:
224 if(len(vis) >0): # (note: in case len is 1, we only copy, essentially)
225 casalog.post('copying '+vis[0]+' to '+theconcatvis , 'INFO')
226 shutil.copytree(vis[0], theconcatvis)
227 # note that the resulting copy is writable even if the original was read-only
228 vis.pop(0)
229 # don't need to pop visweightscale here!
232 if not copypointing: # remove the rows from the POINTING table of the first MS
233 casalog.post('*** copypointing==False: resulting MS will have empty POINTING table.', 'INFO')
234 tmptabname = 'TMPPOINTING'+str(time.time())
235 shutil.rmtree(tmptabname, ignore_errors=True)
236 shutil.move(theconcatvis+'/POINTING', tmptabname)
237 t.open(tmptabname)
238 if(t.nrows()>0):
239 ttab = t.copy(newtablename=theconcatvis+'/POINTING', deep=False, valuecopy=True, norows=True)
240 ttab.close()
241 t.close()
242 shutil.rmtree(tmptabname, ignore_errors=True)
243 else: # the POINTING table is already empty
244 casalog.post('*** Input POINTING table was already empty.', 'INFO')
245 shutil.move(tmptabname, theconcatvis+'/POINTING')
246 t.close()
249 # handle the ephemeris concatenation
250 if not forcesingleephemfield=='':
251 from .concatephem import findephems, concatephem
253 if type(forcesingleephemfield)==str or type(forcesingleephemfield)==int:
254 forcesingleephemfield = [forcesingleephemfield]
255 if not type(forcesingleephemfield) == list:
256 raise RuntimeError('Type of parameter forcesingleephemfield must be str, int, or list')
258 themss = [theconcatvis]
259 for x in vis:
260 themss.append(x)
262 for ephemfield in forcesingleephemfield:
263 if not type(ephemfield)==int:
264 ephemfield = str(ephemfield)
265 casalog.post('*** Forcing single ephemeris for field '+str(ephemfield), 'INFO')
266 thetabs = findephems(themss, ephemfield)
267 if thetabs != [] and not ('' in thetabs):
268 tmptab = os.path.basename(thetabs[0])+'.concattmp'
269 targettab = theconcatvis+'/FIELD/'+os.path.basename(thetabs[0])
270 if not os.path.exists(targettab):
271 raise RuntimeError('Internal ERROR: ephemeris '+targettab+' does not exist')
272 concatephem(thetabs, tmptab)
273 if os.path.exists(tmptab):
274 os.system('rm -rf '+targettab)
275 os.system('mv '+tmptab+' '+targettab)
276 else:
277 casalog.post('ERROR while forcing single ephemeris for field '+str(ephemfield), 'SEVERE')
278 raise RuntimeError('Concatenation of ephemerides for field '+str(ephemfield)+' failed.')
279 else:
280 casalog.post('ERROR while forcing single ephemeris for field '+str(ephemfield), 'SEVERE')
281 raise RuntimeError('Cannot find ephemerides for field '+str(ephemfield)+' in all input MSs.')
284 # Determine if scratch columns should be considered at all
285 # by checking if any of the MSs has them.
287 considerscrcols = False
288 considercorr = False
289 considermodel = False
290 needscrcols = []
291 needmodel = []
292 needcorr = []
293 if ((type(theconcatvis)==str) and (os.path.exists(theconcatvis))):
295 # check if all scratch columns are present
296 t.open(theconcatvis)
297 if (t.colnames().count('MODEL_DATA')==1):
298 considermodel = True
299 if(t.colnames().count('CORRECTED_DATA')==1):
300 considercorr = True
302 needscrcols.append(t.colnames().count('CORRECTED_DATA')==0
303 or t.colnames().count('MODEL_DATA')==0)
304 needmodel.append(t.colnames().count('MODEL_DATA')==0)
305 needcorr.append(t.colnames().count('CORRECTED_DATA')==0)
306 t.close()
307 else:
308 raise ValueError('Visibility data set '+theconcatvis+' not found - please verify the name')
311 for elvis in vis : ###Oh no Elvis does not exist Mr Bill
312 if(not os.path.exists(elvis)):
313 raise ValueError('Visibility data set '+elvis+' not found - please verify the name')
315 # check if all scratch columns are present
316 t.open(elvis)
317 if (t.colnames().count('MODEL_DATA')==1):
318 considermodel = True
319 if(t.colnames().count('CORRECTED_DATA')==1):
320 considercorr = True
322 needscrcols.append(t.colnames().count('CORRECTED_DATA')==0
323 or t.colnames().count('MODEL_DATA')==0)
324 needmodel.append(t.colnames().count('MODEL_DATA')==0)
325 needcorr.append(t.colnames().count('CORRECTED_DATA')==0)
326 t.close()
328 considerscrcols = (considercorr or considermodel) # there are scratch columns
331 # start actual work, file existence has already been checked
332 i = 0
333 if(considerscrcols and needscrcols[i]):
334 # create scratch cols
335 casalog.post('creating scratch columns in '+theconcatvis , 'INFO')
336 _cb.open(theconcatvis,
337 addcorr=(considercorr and needcorr[i]),
338 addmodel=(considermodel and needmodel[i])) # calibrator-open creates scratch columns
339 _cb.close()
341 # scale the weights and sigma of the first MS in the chain
342 if doweightscale:
343 wscale = visweightscale[i]
344 if(wscale==1.):
345 casalog.post('Will leave the weights for this MS unchanged.', 'INFO')
346 else:
347 casalog.post('Scaling weights for first MS by factor '+str(wscale), 'INFO')
348 t.open(theconcatvis, nomodify=False)
349 for colname in [ 'WEIGHT', 'WEIGHT_SPECTRUM']:
350 if (colname in t.colnames()) and (t.iscelldefined(colname,0)):
351 for j in range(0,t.nrows()):
352 a = t.getcell(colname, j)
353 a *= wscale
354 t.putcell(colname, j, a)
355 for colname in ['SIGMA']:
356 if (wscale > 0. and colname in t.colnames()) and (t.iscelldefined(colname,0)):
357 sscale = 1./sqrt(wscale)
358 for j in range(0,t.nrows()):
359 a = t.getcell(colname, j)
360 a *= sscale
361 t.putcell(colname, j, a)
362 t.close()
364 # determine handling switch value
365 handlingswitch = 0
366 if not copypointing:
367 handlingswitch = 2
369 m.open(theconcatvis,nomodify=False)
370 mmsmembers = [theconcatvis]
372 for elvis in vis :
373 i = i + 1
374 destms = ""
375 casalog.post('concatenating '+elvis+' into '+theconcatvis , 'INFO')
377 wscale = 1.
378 if doweightscale:
379 wscale = visweightscale[i]
380 if(wscale==1.):
381 casalog.post('Will leave the weights for this MS unchanged.', 'INFO')
382 else:
383 casalog.post('Will scale weights for this MS by factor '+str(wscale) , 'INFO')
385 if(considerscrcols and needscrcols[i]):
386 if(ParallelTaskHelper.isParallelMS(elvis)):
387 raise RuntimeError('Cannot create scratch columns in a multi-MS. Use virtualconcat.')
388 else:
389 # create scratch cols
390 casalog.post('creating scratch columns for '+elvis+' (original MS unchanged)', 'INFO')
391 tempname = elvis+'_with_scrcols'
392 shutil.rmtree(tempname, ignore_errors=True)
393 shutil.copytree(elvis, tempname)
394 _cb.open(tempname,
395 addcorr=(considercorr and needcorr[i]),
396 addmodel=(considermodel and needmodel[i])) # calibrator-open creates scratch columns
397 _cb.close()
398 # concatenate copy instead of original file
399 m.concatenate(msfile=tempname,freqtol=freqtol,dirtol=dirtol,respectname=respectname,
400 weightscale=wscale,handling=handlingswitch,
401 destmsfile=destms)
402 shutil.rmtree(tempname, ignore_errors=True)
403 else:
404 m.concatenate(msfile=elvis,freqtol=freqtol,dirtol=dirtol,respectname=respectname,
405 weightscale=wscale,handling=handlingswitch,
406 destmsfile=destms)
408 if timesort:
409 casalog.post('Sorting main table by TIME ...', 'INFO')
410 m.timesort()
412 # Write history to output MS, not the input ms.
413 try:
414 param_names = concat.__code__.co_varnames[:concat.__code__.co_argcount]
415 vars = locals( )
416 param_vals = [vars[p] for p in param_names]
417 write_history(mstool(), concatvis, 'concat', param_names,
418 param_vals, casalog)
420 except Exception as instance:
421 casalog.post("*** Error \'%s\' updating HISTORY" % (instance),
422 'WARN')
425 finally:
426 m.close()