Coverage for /wheeldirectory/casa-6.7.0-12-py3.10.el8/lib/py/lib/python3.10/site-packages/casatasks/private/concatephem.py: 58%

1import os 

2import glob 

3 

4from casatools import table as tbtool 

5from casatasks import casalog 

6 

7 

8# Concatenation of ephemeris tables 

9# 

10# Example: 

11# import recipes.ephemerides.concatephem as cce 

12# cce.concatephem(ephems=['Neptune1.tab', 'Neptune2.tab', 'Neptune3.tab'], 

13# outputephem='NeptuneAll.tab') 

14 

15 

16def concatephem(ephems=[], outputephem=''): 

17 """ 

18 concatephem 

19 

20 Concatenate the given ephemeris tables into a single output table 

21 filling gaps with dummy rows and removing overlap rows. 

22 Before concatenation, test that basic conditions are fulfilled: 

23 same time grid, list of input ephemerides is ordered in time. 

24 

25 ephems - List of the ephemeris tables to be concatenated 

26 default:  

27 outputephem - Name of the output ephemeris to be created from the concatenation 

28 If empty, concatephem will only perform a dryrun. 

29 default: '' 

30  

31 """ 

32 mytb = tbtool() 

33 starts = [] 

34 ends = [] 

35 stepsizes = [] 

36 hasoverlap_with_previous = [] 

37 gap_to_previous = [] # negative entry means overlap  

38 shouldconcat_with_previous = [] 

39 canconcat_with_previous = [] 

40 

41 stepsize_rel_tolerance = 1E-6 # relative difference in time step size must be smaller than this 

42 stepsize_abs_tolerance_d = 0.1/86400. # absolute difference in time step size (days) must be smaller than this 

43 

44 for myephem in ephems: 

45 mytb.open(myephem) 

46 mynrows = mytb.nrows() 

47 if mynrows==0: 

48 mytb.close() 

49 casalog.post('Ephemeris '+myephem+' has no rows.', 'SEVERE') 

50 return 

51 mystart = mytb.getcell('MJD',0) 

52 myend = mytb.getcell('MJD', mynrows-1) 

53 if mynrows<2: 

54 mystepsize = 0 

55 casalog.post('Ephemeris '+myephem+' has only one row.', 'WARN') 

56 else: 

57 mystepsize = mytb.getcell('MJD',1) - mystart 

58 

59 mytb.close() 

60 

61 starts.append(mystart) 

62 ends.append(myend) 

63 stepsizes.append(mystepsize) 

64 

65 if os.path.exists(outputephem): 

66 casalog.post('Output ephemeris table '+outputephem+' exists already.', 'SEVERE') 

67 return 

68 

69 casalog.post('Ephem no., startMJD, endMJD, step size (d)', 'INFO') 

70 for i in range(0,len(starts)): 

71 casalog.post(str(i)+', '+str(starts[i])+', '+str(ends[i])+', '+str(stepsizes[i]), 'INFO') 

72 

73 backstep = 0 

74 for i in range(0,len(starts)): 

75 shouldconcat_with_previous.append(False) 

76 canconcat_with_previous.append(False) 

77 hasoverlap_with_previous.append(False) 

78 gap_to_previous.append(0) 

79 if i>0: 

80 if (abs(stepsizes[i] - stepsizes[i-1-backstep])/stepsizes[i] < stepsize_rel_tolerance) \ 

81 and abs(stepsizes[i] - stepsizes[i-1-backstep]) < stepsize_abs_tolerance_d: 

82 casalog.post( 'Ephemerides '+str(i-1-backstep)+' and '+str(i)+' have same step size.', 'INFO') 

83 if starts[i-1-backstep] <= starts[i]: 

84 casalog.post( 'Ephemeris '+str(i-1-backstep)+' begins before '+str(i), 'INFO') 

85 if ends[i-1-backstep] < ends[i]: 

86 casalog.post( 'Ephemeris '+str(i-1-backstep)+' ends before '+str(i), 'INFO') 

87 shouldconcat_with_previous[i] = True 

88 numsteps_to_add = (starts[i]-ends[i-1-backstep])/stepsizes[i] - 1 

89 gap_to_previous[i] = int(round(numsteps_to_add)) 

90 if abs(round(numsteps_to_add) - numsteps_to_add)<1E-4: 

91 casalog.post( 'Gap between ephemerides '+str(i-1-backstep)+' and '+str(i)+' is '+str(gap_to_previous[i])+' steps', 'INFO') 

92 canconcat_with_previous[i] = True 

93 backstep = 0 

94 if numsteps_to_add<0.: 

95 casalog.post( 'Ephemerides '+str(i-1-backstep)+' and '+str(i)+' have overlap of '+str(-gap_to_previous[i])+' steps', 'INFO') 

96 hasoverlap_with_previous[i]=True 

97 else: 

98 casalog.post( 'Gap between ephemerides '+str(i-1-backstep)+' and '+str(i)+' is not an integer number of steps:', 'WARN') 

99 casalog.post( str(round(numsteps_to_add) - numsteps_to_add), 'INFO') 

100 canconcat_with_previous[i] = False 

101 backstep += 1 

102 else: 

103 casalog.post( 'Ephemeris '+str(i-1-backstep)+' does not end before '+str(i), 'INFO') 

104 shouldconcat_with_previous[i] = False 

105 canconcat_with_previous[i] = False 

106 backstep += 1 

107 else: 

108 casalog.post( 'Ephemeris '+str(i-1-backstep)+' does not begin before or at the same time as '+str(i), 'INFO') 

109 casalog.post('Ephemerides are in wrong order.', 'SEVERE') 

110 return 

111 else: 

112 casalog.post( 'Ephemerides '+str(i-1-backstep)+' and '+str(i)+' do not have same step size.', 'INFO') 

113 casalog.post('Ephemerides have inhomogeneous steps sizes.', 'SEVERE') 

114 return 

115 #end if 

116 

117 

118 if outputephem=='' or len(starts)<2: 

119 return 

120 else: 

121 casalog.post( 'Creating output ephemeris ...', 'INFO') 

122 

123 os.system('cp -R '+ephems[0]+' '+outputephem) 

124 

125 mytb2 = tbtool() 

126 

127 try: 

128 mytb.open(outputephem, nomodify=False) 

129 except: 

130 casalog.post('Error opening table '+outputepehem+' for writing.', 'SEVERE') 

131 return 

132 

133 

134 for i in range(1,len(starts)): 

135 

136 if shouldconcat_with_previous[i] and canconcat_with_previous[i]: 

137 

138 mynrows = mytb.nrows() 

139 

140 mytb2.open(ephems[i]) 

141 mynrows2 = mytb2.nrows() 

142 startrow = 0 

143 if(hasoverlap_with_previous[i]): 

144 startrow = -gap_to_previous[i] 

145 elif(gap_to_previous[i]>0): # fill gap with dummy rows 

146 mytb.addrows(gap_to_previous[i]) 

147 startmjd = mytb.getcell('MJD', mynrows-1) 

148 for j in range(mynrows, mynrows+gap_to_previous[i]): 

149 newmjd = startmjd+stepsizes[i]*(j-mynrows+1) 

150 mytb.putcell('MJD', j, newmjd) 

151 mytb.putcell('RA', j, 0.) 

152 mytb.putcell('DEC', j, 0.) 

153 mytb.putcell('Rho', j, 0.) 

154 mytb.putcell('RadVel', j, 0.) 

155 mytb.putcell('diskLong', j, 0.) 

156 mytb.putcell('diskLat', j, 0.) 

157 

158 casalog.post( str(gap_to_previous[i])+' dummy rows appended to fill gap', 'INFO') 

159 

160 mynrows = mytb.nrows() 

161 

162 mytb.addrows(mynrows2-startrow) 

163 for j in range(mynrows, mynrows+mynrows2-startrow): 

164 fromrow = j - mynrows + startrow 

165 mytb.putcell('MJD', j, mytb2.getcell('MJD', fromrow)) 

166 mytb.putcell('RA', j, mytb2.getcell('RA', fromrow)) 

167 mytb.putcell('DEC', j, mytb2.getcell('DEC', fromrow)) 

168 mytb.putcell('Rho', j, mytb2.getcell('Rho', fromrow)) 

169 mytb.putcell('RadVel', j, mytb2.getcell('RadVel', fromrow)) 

170 mytb.putcell('diskLong', j, mytb2.getcell('diskLong', fromrow)) 

171 mytb.putcell('diskLat', j, mytb2.getcell('diskLat', fromrow)) 

172 

173 casalog.post( str(mynrows2-startrow)+' rows copied over', 'INFO') 

174 

175 else: 

176 

177 casalog.post( 'Ephemeris '+str(i)+' skipped', 'INFO') 

178 

179 #endif 

180 #endfor 

181 #endif 

182 

183 return 

184 

185 

186def findephems(vis=[], field=''): 

187 

188 """ 

189 findephems 

190 

191 Search the MSs given in the list "vis" for ephemerides for 

192 a given field and return the list of paths to the ephemeris tables 

193 in the same order as vis. 

194 

195 vis - list of the MSs to search for ephemerides 

196 default: [] 

197 

198 field - field for which to seach ephemerides 

199 default: 

200 """ 

201 

202 if vis==[] or field=='': 

203 return [] 

204 

205 if type(vis) == str: 

206 vis = [vis] 

207 

208 mytb=tbtool() 

209 

210 rval = [] # list of ephemeris tables to be returned 

211 

212 for visname in vis: 

213 if not os.path.exists(visname+'/FIELD'): 

214 casalog.post('MS '+visname+' does not exist.', 'SEVERE') 

215 return [] 

216 

217 mytb.open(visname+'/FIELD') 

218 fnames = mytb.getcol('NAME') 

219 

220 thefield = field 

221 if type(thefield) == int: 

222 thefield = str(fnames[thefield]) 

223 if type(thefield) != str: 

224 casalog.post('Parameter field must be str or int.', 'SEVERE') 

225 mytb.close() 

226 return [] 

227 

228 if thefield in fnames: 

229 colnames = mytb.colnames() 

230 if not 'EPHEMERIS_ID' in colnames: 

231 casalog.post('MS '+visname+' has no ephemerides attached.', 'WARN') 

232 rval.append('') 

233 mytb.close() 

234 continue 

235 

236 ephids = mytb.getcol('EPHEMERIS_ID') 

237 mytb.close() 

238 theephid = ephids[list(fnames).index(thefield)] 

239 

240 thetabs = glob.glob(visname+'/FIELD/EPHEM'+str(theephid)+'_*') 

241 

242 if len(thetabs)==0: 

243 casalog.post('MS '+visname+' has no ephemerides for field '+thefield+' attached.', 'WARN') 

244 rval.append('') 

245 continue 

246 

247 if len(thetabs)>1: 

248 casalog.post('MS '+visname+' has more than one ephemeris for field '+thefield+' attached.', 'SEVERE') 

249 return[] 

250 

251 rval.append(thetabs[0]) 

252 

253 else: 

254 casalog.post('MS '+visname+' has no field '+thefield, 'WARN') 

255 rval.append('') 

256 continue 

257 #endfor 

258 

259 return rval 

260 

261 

262def concatreplaceephem(vis=[], field='', commontime=False): 

263 

264 """ 

265 Search the MSs given in the list "vis" for ephemerides for 

266 a given field, concatenate the ephemerides, and replace 

267 all of the original ephemerides with the concatenated one. 

268 

269 vis - list of the MSs to search for ephemerides 

270 default: [] 

271 

272 field - field for which to seach ephemerides 

273 default: 

274 

275 commontime - set the FIELD table reference time in all input MSs to a common time. 

276 If True, the common time is taken from the time of the first MS. 

277 If a float value is provided, it is interpreted as a time in seconds in the native 

278 CASA time representation. 

279 default: False - do not modify the FIELD table TIME column 

280 

281 """ 

282 

283 ephemfield = field 

284 thetabs = findephems(vis, ephemfield) 

285 

286 if not (type(commontime)==bool) and not (type(commontime)==float): 

287 raise Exception('ERROR: parameter commontime must be bool or float') 

288 

289 if thetabs != [] and not ('' in thetabs): 

290 tmptab = os.path.basename(thetabs[0])+'.concattmp' 

291 concatephem(thetabs, tmptab) 

292 if os.path.exists(tmptab): 

293 for targettab in thetabs: 

294 if not os.path.exists(targettab): 

295 raise Exception('Internal ERROR: ephemeris '+targettab+' does not exist') 

296 os.system('rm -rf '+targettab) 

297 os.system('cp -R '+tmptab+' '+targettab) 

298 else: 

299 casalog.post('ERROR while concatenating ephemerides for field '+str(ephemfield), 'SEVERE') 

300 raise Exception('Concatenation of ephemerides for field '+str(ephemfield)+' failed.') 

301 

302 os.system('rm -rf '+tmptab) 

303 

304 if type(commontime)==float or commontime==True: 

305 casalog.post('Will set FIELD table reference time for field '+str(ephemfield)+' to common value:', 'INFO') 

306 thecommontime = commontime 

307 if type(commontime)==float: 

308 casalog.post(' '+str(thecommontime), 'INFO') 

309 mytb = tbtool() 

310 for thevis in vis: 

311 mytb.open(thevis+'/FIELD', nomodify=False) 

312 thenames = mytb.getcol('NAME') 

313 thetimes = mytb.getcol('TIME') 

314 

315 for i in range(0,len(thenames)): 

316 if (thenames[i]==ephemfield): 

317 if thecommontime==True and thevis==vis[0]: 

318 thecommontime = thetimes[i] # memorize the common time 

319 casalog.post(' '+str(thecommontime), 'INFO') 

320 else: 

321 thetimes[i] = thecommontime 

322 #end for 

323 mytb.putcol('TIME', thetimes) 

324 #end for 

325 mytb.close() 

326 

327 else: 

328 casalog.post('ERROR while searching for ephemerides for field '+str(ephemfield), 'SEVERE') 

329 raise Exception('Cannot find ephemerides for field '+str(ephemfield)+' in all input MSs.') 

330 

331 casalog.post('All ephemerides for field '+str(ephemfield)+' replaced by concatenated ephemeris.', 'INFO') 

332 

333 return True