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

274 statements  

« prev     ^ index     » next       coverage.py v7.6.4, created at 2024-11-01 07:19 +0000

1import collections 

2import contextlib 

3import itertools 

4import os 

5import shutil 

6 

7import numpy as np 

8 

9from casatasks import casalog 

10from casatasks.private import sdutil, simutil 

11from casatools import ms as mstool 

12from casatools import msmetadata, quanta, singledishms 

13 

14ut = simutil.simutil() 

15qa = quanta() 

16sdms = singledishms() 

17 

18 

19class ATMParameterConfigurator(collections.abc.Iterator): 

20 def __init__(self, key, value, do_config=True): 

21 data = [(key, value)] if do_config else [] 

22 self._iter = iter(data) 

23 

24 def __next__(self): 

25 return next(self._iter) 

26 

27 

28@contextlib.contextmanager 

29def open_msmd(path): 

30 msmd = msmetadata() 

31 msmd.open(path) 

32 try: 

33 yield msmd 

34 finally: 

35 msmd.close() 

36 

37 

38def _ms_remove(path): 

39 if (os.path.exists(path)): 

40 if (os.path.isdir(path)): 

41 shutil.rmtree(path) 

42 else: 

43 os.remove(path) 

44 

45 

46def get_default_params(): 

47 # Default constant: taken from atmcor_20200807.py (CSV-3320) 

48 atmtype = 2 # atmType parameter for at (1: tropical, 2: mid lat summer, 3: mid lat winter, etc) 

49 maxalt = 120 # maxAltitude parameter for at (km) 

50 lapserate = -5.6 # dTem_dh parameter for at (lapse rate; K/km) 

51 scaleht = 2.0 # h0 parameter for at (water scale height; km) 

52 

53 # convolve dTa* spectra with [0.25, 0.5, 0.25] to mimic Hanning spectral response; 

54 dosmooth = False 

55 # set to True if spectral averaging was not employed for the spw 

56 dp = 10.0 # initATMProfile DEFAULT ### 

57 dpm = 1.2 # initATMProfile DEFAULT ### 

58 return locals() 

59 

60 

61def parse_gainfactor(gainfactor): 

62 """Parse gainfactor parameter. 

63 

64 Parse gainfactor parameter. 

65 

66 Args: 

67 gainfactor (float, dict, str): gain factor. 

68 if float value is given, it applies to all spws. 

69 if dict is given, spw id and corresponding factor 

70 should be provided as key-value pair. 

71 if str is given, it should be the name of caltable. 

72 factors are derived as inverse-square of values 

73 stored in the caltable. 

74 

75 Raises: 

76 FileNotFoundError: specified caltable does not exist. 

77 

78 Returns: 

79 dictionary whose keys are spw id in string while values 

80 are the factors to be applied to each spw. 

81 dictionary is defined as collections.defaultdict that 

82 returns 1.0 as a default value. 

83 """ 

84 gaindict = collections.defaultdict(lambda: 1.0) 

85 if isinstance(gainfactor, dict): 

86 # make sure keys are str 

87 d = dict((str(k), v) for k, v in gainfactor.items()) 

88 gaindict.update(d) 

89 elif isinstance(gainfactor, str): 

90 # should be the name of caltable 

91 if not os.path.exists(gainfactor): 

92 raise FileNotFoundError('"{}" should exist.'.format(gainfactor)) 

93 with sdutil.table_manager(gainfactor) as tb: 

94 if 'FPARAM' in tb.colnames(): 

95 col = 'FPARAM' 

96 elif 'CPARAM' in tb.colnames(): 

97 col = 'CPARAM' 

98 else: 

99 raise RuntimeError('{} is not a caltable'.format(gainfactor)) 

100 spw_list = set(tb.getcol('SPECTRAL_WINDOW_ID')) 

101 for spw in spw_list: 

102 tsel = tb.query('SPECTRAL_WINDOW_ID=={}'.format(spw)) 

103 try: 

104 v = tsel.getcol(col).real 

105 finally: 

106 tsel.close() 

107 factor = np.mean(1 / np.square(v)) 

108 gaindict[str(spw)] = factor 

109 else: 

110 # should be float 

111 v = float(gainfactor) 

112 gaindict = collections.defaultdict(lambda: v) 

113 return gaindict 

114 

115 

116def gaindict2list(msname, gaindict): 

117 with sdutil.table_manager(os.path.join(msname, 'SPECTRAL_WINDOW')) as tb: 

118 nspw = tb.nrows() 

119 

120 gainlist = np.ones(nspw, dtype=float) 

121 if isinstance(gaindict, collections.defaultdict) and len(gaindict.keys()) == 0: 

122 gainlist[:] = gaindict[0] 

123 else: 

124 for k, v in gaindict.items(): 

125 spw = int(k) 

126 if 0 <= spw and spw < nspw: 

127 gainlist[spw] = v 

128 

129 return gainlist 

130 

131 

132def get_all_spws_from_main(msname): 

133 """Extract all spw ids from MAIN table. 

134 

135 Extract all spectral window ids that have any 

136 associated data in MS MAIN table. 

137 

138 Args: 

139 msname (str): name of MS 

140 

141 Returns: 

142 list: list of available spectral window ids 

143 """ 

144 with sdutil.table_manager(msname) as tb: 

145 ddids = np.unique(tb.getcol('DATA_DESC_ID')) 

146 with open_msmd(msname) as msmd: 

147 spws_all = [msmd.spwfordatadesc(ddid) for ddid in ddids] 

148 return spws_all 

149 

150 

151def get_selected_spws(msname, spw): 

152 """ 

153 Get selected spectral window ids. 

154 

155 Args: 

156 msname (str): name of MS 

157 spw (str): spectral window selection 

158 

159 Raises: 

160 TypeError: spw is not string 

161 

162 Returns: 

163 list: list of selected spectral window ids 

164 """ 

165 if not isinstance(spw, str): 

166 raise TypeError('spw selection must be string') 

167 elif len(spw) == 0: 

168 # '' indicates all spws, which is equivalent to '*' 

169 spwsel = '*' 

170 else: 

171 spwsel = spw 

172 ms = mstool() 

173 sel = ms.msseltoindex(msname, spw=spwsel) 

174 return sel['spw'] 

175 

176 

177def parse_spw(msname, spw=''): 

178 """Parse spw selection into list of spw ids. 

179 

180 Parse spw selection into list of spw ids that have 

181 associated data in the MAIN table of given MS. 

182 

183 Args: 

184 msname (str): name of MS 

185 spw (str): spw selection 

186 

187 Raises: 

188 TypeError: spw selection is not str 

189 RuntimeError: spw selection cause empty result 

190 

191 Returns: 

192 list: list of selected spw ids 

193 """ 

194 spws_all = get_all_spws_from_main(msname) 

195 spws_sel = get_selected_spws(msname, spw) 

196 spws = set(spws_all).intersection(set(spws_sel)) 

197 return list(spws) 

198 

199 

200def get_mount_off_source_commands(msname): 

201 """Return list of flag commands whose reason is "Mount_is_off_source". 

202 

203 Args: 

204 msname (str): name of MS 

205 

206 Returns: 

207 np.ndarray: list of flag commands 

208 """ 

209 with sdutil.table_manager(os.path.join(msname, 'FLAG_CMD')) as tb: 

210 if tb.nrows() > 0: 

211 tsel = tb.query('REASON=="Mount_is_off_source"') 

212 try: 

213 commands = tsel.getcol('COMMAND') 

214 finally: 

215 tsel.close() 

216 else: 

217 commands = [] 

218 return commands 

219 

220 

221def get_antenna_name(antenna_selection): 

222 """Extract antenna name from the antenna selection string. 

223 

224 Here, antenna_selection is assumed to be a string 

225 in the form '<ANTENNA_NAME>&&*'. 

226 

227 Args: 

228 antenna_selection (str): antenna selection string 

229 

230 Returns: 

231 str: antenna name 

232 """ 

233 return antenna_selection.split('=')[1].strip("'&*") 

234 

235 

236def get_time_delta(time_range): 

237 """Convert time range string into time duration in sec. 

238 

239 Here, time_range is assumed to be a string in the form 

240 'YYYY/MM/DD/hh:mm:ss~YYYY/MM/DD/hh:mm:ss' 

241 

242 Args: 

243 time_range (str): time range string 

244 

245 Returns: 

246 float: time duration in sec 

247 """ 

248 timestrs = time_range.split('=')[1].strip("'").split('~') 

249 timequanta = [qa.quantity(t) for t in timestrs] 

250 timedelta = qa.convert(qa.sub(timequanta[1], timequanta[0]), 's')['value'] 

251 return abs(timedelta) 

252 

253 

254def cmd_to_ant_and_time(cmd): 

255 """Extract antenna name and time duration from the flag command. 

256 

257 Args: 

258 cmd (str): flag command 

259 

260 Returns: 

261 tuple: antenna name and time duration 

262 """ 

263 sels = cmd.split() 

264 asel = list(filter(lambda x: x.startswith('antenna'), sels))[0] 

265 tsel = list(filter(lambda x: x.startswith('time'), sels))[0] 

266 

267 antenna_name = get_antenna_name(asel) 

268 time_delta = get_time_delta(tsel) 

269 

270 return antenna_name, time_delta 

271 

272 

273def inspect_flag_cmd(msname): 

274 """Inspect FLAG_CMD table. 

275 

276 Search flag commands whose reason is Mount_is_off_source and 

277 extract antenna name and time duration from the commands. 

278 

279 Args: 

280 msname (str): name of MS 

281 

282 Returns: 

283 tuple: two dictionaries containing the inspection result. 

284 The first one is number of command counts per antenna 

285 while the second one is total duration flagged by the 

286 commands per antenna. 

287 """ 

288 commands = get_mount_off_source_commands(msname) 

289 

290 cmd_counts = collections.defaultdict(lambda: 0) 

291 time_counts = collections.defaultdict(lambda: 0) 

292 

293 for cmd in commands: 

294 ant, dt = cmd_to_ant_and_time(cmd) 

295 cmd_counts[ant] += 1 

296 time_counts[ant] += dt 

297 

298 return cmd_counts, time_counts 

299 

300 

301# Argument parameter handling 

302def parse_atm_params(user_param, user_default, task_default, default_unit=''): 

303 """Parse ATM parameters. 

304 

305 Args: 

306 user_param (str,int,float): User input. 

307 user_default (str,int,float): User default. 

308 task_default (str,int,float): Task default. 

309 default_unit (str): Default unit. 

310 

311 Raises: 

312 ValueError: user_param is invalid. 

313 

314 Returns: 

315 Tuple: Two-tuple, resulting value as quantity and boolean 

316 value indicating if the value is equal to user_default. 

317 """ 

318 is_customized = user_param != user_default and user_param is not None 

319 

320 try: 

321 if qa.isquantity(task_default): 

322 task_default_quanta = qa.quantity(task_default) 

323 else: 

324 task_default_quanta = qa.quantity(task_default, default_unit) 

325 except Exception as e: 

326 casalog.post('INTERNAL ERROR: {}'.format(e), priority='SEVERE') 

327 raise 

328 

329 if not is_customized: 

330 param = task_default_quanta['value'] 

331 else: 

332 user_param_quanta = qa.quantity(user_param) 

333 if user_param_quanta['unit'] == '': 

334 user_param_quanta = qa.quantity( 

335 user_param_quanta['value'], 

336 default_unit 

337 ) 

338 else: 

339 user_param_quanta = qa.convert( 

340 user_param_quanta, 

341 default_unit 

342 ) 

343 is_compatible = qa.compare(user_param_quanta, task_default_quanta) 

344 if is_compatible: 

345 param = user_param_quanta['value'] 

346 else: 

347 raise ValueError('User input "{}" should have the unit compatible with "{}"'.format( 

348 user_param, 

349 default_unit 

350 )) 

351 

352 return param, is_customized 

353 

354 

355def parse_atm_list_params(user_param, user_default='', task_default=[], default_unit=''): 

356 """Parse ATM parameters. 

357 

358 Args: 

359 user_param (str,list): User input. 

360 user_default (str): User default. 

361 task_default (list): Task default. 

362 default_unit (str): Unit for output values. 

363 

364 Raises: 

365 ValueError: user_param is invalid. 

366 

367 Returns: 

368 Tuple: Two-tuple, resulting value as quantity and boolean 

369 value indicating if the value is equal to user_default. 

370 """ 

371 is_customized = user_param != user_default and user_param is not None 

372 

373 if not is_customized: 

374 return task_default, is_customized 

375 

376 if isinstance(user_param, (list, np.ndarray)): 

377 try: 

378 param = [parse_atm_params(p, user_default, 0, default_unit=default_unit)[0] 

379 for p in user_param] 

380 param = [qa.convert(p, default_unit)['value'] for p in param] 

381 except Exception as e: 

382 casalog.post('ERROR during handling list input: {}'.format(e)) 

383 raise ValueError('list input "{}" is invalid.'.format(user_param)) 

384 return param, is_customized 

385 elif isinstance(user_param, str): 

386 try: 

387 split_param = user_param.split(',') 

388 param, _ = parse_atm_list_params(split_param, user_default, task_default, default_unit) 

389 except Exception as e: 

390 casalog.post('ERROR during handling comma-separated str input: {}'.format(e)) 

391 raise ValueError('str input "{}" is invalid.'.format(user_param)) 

392 return param, is_customized 

393 else: 

394 raise ValueError('user_param for parse_atm_list_params should be either list or str.') 

395 

396 

397def get_default_antenna(msname): 

398 """Determine default antenna id based on the FLAG_CMD table. 

399 

400 Procedure is as follows. 

401 

402 (1) extract flag commands whose reason is "Mount_is_off_source". 

403 (2) compile the commands into a number of commands and flagged 

404 time durations for each antenna. 

405 (3) select antenna with the shortest flagged duration. 

406 (4) if multiple antennas match in (3), select antenna with 

407 the least number of commands among them. 

408 (5) if multiple antennas match in (4), select the first 

409 antenna among them. 

410 

411 Args: 

412 msname (str): name of MS 

413 

414 Raises: 

415 Exception: no antenna was found in the MAIN table 

416 

417 Returns: 

418 int: default antenna id 

419 """ 

420 # get list of antenna Ids from MAIN table 

421 with sdutil.table_manager(msname) as tb: 

422 ant_list = np.unique(tb.getcol('ANTENNA1')) 

423 

424 # No Available antenna 

425 if len(ant_list) == 0: 

426 raise Exception("No Antenna was found.") 

427 

428 # get antenna names list by antenna Id 

429 with open_msmd(msname) as msmd: 

430 ant_name = [msmd.antennanames(i)[0] for i in ant_list] 

431 

432 # dictionary to map antenna name to antenna Id 

433 ant_dict = dict((k, v) for k, v in zip(ant_name, ant_list)) 

434 

435 # determine default antenna id 

436 cmd_counts, flagged_durations = inspect_flag_cmd(msname) 

437 

438 if len(cmd_counts) == 0: 

439 # No flag command exists. All the antennas should be healthy 

440 # so just pick up the first antenna. 

441 default_id = ant_list[0] 

442 default_name = ant_name[0] 

443 else: 

444 flagged_durations_filtered = dict((k, flagged_durations[k]) for k in ant_dict.keys()) 

445 min_duration = min(flagged_durations_filtered.values()) 

446 candidate_antennas = [k for k, v in flagged_durations_filtered.items() if v == min_duration] 

447 

448 if len(candidate_antennas) == 1: 

449 default_name = candidate_antennas[0] 

450 default_id = ant_dict[default_name] 

451 else: 

452 _counts = [cmd_counts[a] for a in candidate_antennas] 

453 min_count = min(_counts) 

454 candidate_antennas2 = [a for i, a in enumerate(candidate_antennas) 

455 if _counts[i] == min_count] 

456 default_name = candidate_antennas2[0] 

457 default_id = ant_dict[default_name] 

458 casalog.post('Select {} (ID {}) as a default antenna'.format(default_name, default_id)) 

459 return default_id 

460 

461 

462def get_default_altitude(msname, antid): 

463 """Get default altitude of the antenna. 

464 

465 decide default value of 'Altitude' for Atm Correction. 

466 This requires to calculate Elevation from Antenna Position Information. 

467 """ 

468 with sdutil.table_manager(os.path.join(msname, 'ANTENNA')) as tb: 

469 # obtain the antenna Position (Earth Center) specified by antid 

470 X, Y, Z = (float(i) for i in tb.getcell('POSITION', antid)) 

471 

472 # xyz2long() -- https://casa.nrao.edu/casadocs/casa-5.6.0/simulation/simutil 

473 # 

474 # When given ITRF Earth-centered (X, Y, Z, using the parameters x, y, and z) 

475 # coordinates [m] for a point, 

476 # this method returns geodetic latitude and longitude [radians] and elevation [m]. 

477 # Elevation is measured relative to the closest point to the (latitude, longitude) 

478 # on the WGS84 (World Geodetic System 1984) reference ellipsoid. 

479 

480 # [0]:longitude, [1]:latitude, [2]:elevation (geodetic elevation) 

481 P = ut.xyz2long(X, Y, Z, 'WGS84') 

482 geodetic_elevation = P[2] 

483 

484 ref = tb.getcolkeyword('POSITION', 'MEASINFO')['Ref'] 

485 

486 casalog.post("Default Altitude") 

487 casalog.post(" - Antenna ID: %d. " % antid) 

488 casalog.post(" - Ref = %s. " % ref) 

489 casalog.post(" - Position: (%s, %s, %s)." % (X, Y, Z)) 

490 casalog.post(" Altitude (geodetic elevation): %f" % geodetic_elevation) 

491 

492 return geodetic_elevation 

493 

494 

495class ATMScalarParameterConfigurator(ATMParameterConfigurator): 

496 def __init__(self, key, user_input, impl_default, default_unit, 

497 api_default='', is_mandatory=True, is_effective=True): 

498 value, is_customized = parse_atm_params( 

499 user_param=user_input, user_default=api_default, 

500 task_default=impl_default, default_unit=default_unit) 

501 do_config = is_mandatory or (is_effective and is_customized) 

502 super().__init__(key=key, value=value, do_config=do_config) 

503 

504 

505class ATMListParameterConfigurator(ATMParameterConfigurator): 

506 def __init__(self, key, user_input, impl_default, default_unit, 

507 api_default='', is_mandatory=True, is_effective=True): 

508 value, is_customized = parse_atm_list_params( 

509 user_param=user_input, user_default=api_default, 

510 task_default=impl_default, default_unit=default_unit) 

511 do_config = is_mandatory or (is_effective and is_customized) 

512 super().__init__(key=key, value=value, do_config=do_config) 

513 

514 

515def get_configuration_for_atmcor(infile, spw, outputspw, gainfactor, user_inputs): 

516 # requested list of output spws and processing spws 

517 # processing spws are the intersection of these 

518 outputspws_param = parse_spw(infile, outputspw) 

519 spws_param = parse_spw(infile, spw) 

520 all_processing_spws = np.asarray(list(set(spws_param).intersection(set(outputspws_param)))) 

521 

522 # generate gain factor dictionary 

523 gaindict = parse_gainfactor(gainfactor) 

524 gainlist = gaindict2list(infile, gaindict) 

525 

526 # default parameter values (from Tsuyoshi's original script) 

527 default_params = get_default_params() 

528 

529 # reference antenna_id to calculate Azimuth/Elevation 

530 reference_antenna = int(get_default_antenna(infile)) 

531 

532 # altitude of reference antenna 

533 default_altitude = get_default_altitude(infile, reference_antenna) 

534 user_altitude = user_inputs['altitude'] if user_inputs['atmdetail'] else '' 

535 

536 parameters = [ 

537 ATMParameterConfigurator(key='processspw', value=all_processing_spws), 

538 ATMParameterConfigurator(key='gainfactor', value=gainlist), 

539 ATMParameterConfigurator(key='refant', value=reference_antenna), 

540 ATMParameterConfigurator(key='atmType', value=user_inputs['atmtype']), 

541 ATMParameterConfigurator(key='maxAltitude', value=float(default_params['maxalt'])), 

542 ATMScalarParameterConfigurator( 

543 key='lapseRate', user_input=user_inputs['dtem_dh'], 

544 impl_default=default_params['lapserate'], default_unit='K/km', 

545 ), 

546 ATMScalarParameterConfigurator( 

547 key='scaleHeight', user_input=user_inputs['h0'], 

548 impl_default=default_params['scaleht'], default_unit='km', 

549 ), 

550 ATMScalarParameterConfigurator( 

551 key='pressureStep', user_input=user_inputs['dp'], 

552 impl_default=default_params['dp'], default_unit='mbar' 

553 ), 

554 ATMScalarParameterConfigurator( 

555 key='pressureStepFactor', user_input=user_inputs['dpm'], 

556 impl_default=default_params['dpm'], default_unit='', api_default=-1 

557 ), 

558 ATMScalarParameterConfigurator( 

559 key='siteAltitude', user_input=user_altitude, 

560 impl_default=default_altitude, default_unit='m', 

561 api_default='', 

562 is_mandatory=True, is_effective=user_inputs['atmdetail'] 

563 ), 

564 ATMScalarParameterConfigurator( 

565 key='pressure', user_input=user_inputs['pressure'], 

566 impl_default=0, default_unit='mbar', 

567 is_mandatory=False, is_effective=user_inputs['atmdetail'] 

568 ), 

569 ATMScalarParameterConfigurator( 

570 key='temperature', user_input=user_inputs['temperature'], 

571 impl_default=0, default_unit='K', 

572 is_mandatory=False, is_effective=user_inputs['atmdetail'] 

573 ), 

574 ATMScalarParameterConfigurator( 

575 key='humidity', user_input=user_inputs['humidity'], 

576 impl_default=0, default_unit='%', 

577 api_default=-1, 

578 is_mandatory=False, is_effective=user_inputs['atmdetail'] 

579 ), 

580 ATMScalarParameterConfigurator( 

581 key='pwv', user_input=user_inputs['pwv'], 

582 impl_default=0, default_unit='mm', 

583 is_mandatory=False, is_effective=user_inputs['atmdetail'] 

584 ), 

585 ATMListParameterConfigurator( 

586 key='layerBoundaries', user_input=user_inputs['layerboundaries'], 

587 impl_default=[], default_unit='m', 

588 is_mandatory=False, is_effective=user_inputs['atmdetail'] 

589 ), 

590 ATMListParameterConfigurator( 

591 key='layerTemperatures', user_input=user_inputs['layertemperature'], 

592 impl_default=[], default_unit='K', 

593 is_mandatory=False, is_effective=user_inputs['atmdetail'] 

594 ) 

595 ] 

596 

597 config = dict(itertools.chain(*parameters)) 

598 

599 # number of threads for OpenMP 

600 # if config['nthreads'] is set to -1, task will decide number of threads automatically 

601 config['nthreads'] = int(os.getenv('OMP_NUM_THREADS', -1)) 

602 

603 return config 

604 

605 

606@sdutil.sdtask_decorator 

607def sdatmcor( 

608 infile=None, datacolumn=None, outfile=None, overwrite=None, 

609 field=None, spw=None, scan=None, antenna=None, 

610 correlation=None, timerange=None, intent=None, 

611 observation=None, feed=None, msselect=None, 

612 outputspw=None, 

613 gainfactor=None, 

614 dtem_dh=None, h0=None, atmtype=None, 

615 atmdetail=None, 

616 altitude=None, temperature=None, pressure=None, humidity=None, pwv=None, 

617 dp=None, dpm=None, 

618 layerboundaries=None, layertemperature=None): 

619 

620 try: 

621 # Input/Output error check and internal set up. 

622 if infile == '': 

623 errmsg = "infile MUST BE specified." 

624 raise Exception(errmsg) 

625 

626 if outfile == '': 

627 errmsg = "outfile MUST BE specified." 

628 raise Exception(errmsg) 

629 

630 # Protection, in case infile == outfile 

631 if infile == outfile: 

632 errmsg = "You are attempting to write the output on your input file." 

633 raise Exception(errmsg) 

634 

635 # File Info 

636 casalog.post("INPUT/OUTPUT") 

637 casalog.post(" Input MS file = %s " % infile) 

638 casalog.post(" Output MS file = %s " % outfile) 

639 

640 # infile Inaccessible 

641 if not os.path.exists(infile): 

642 errmsg = "Specified infile does not exist." 

643 raise Exception(errmsg) 

644 

645 # outfile Protected 

646 if os.path.exists(outfile): 

647 if overwrite: 

648 casalog.post("Overwrite: Overwrite specified. Delete the existing output file.") 

649 _ms_remove(outfile) 

650 else: 

651 errmsg = "Specified outfile already exist." 

652 raise Exception(errmsg) 

653 

654 # Inspect atmtype 

655 atmtype_int = int(atmtype) 

656 if atmtype_int not in (1, 2, 3, 4, 5): 

657 errmsg = "atmtype (=%s) should be any one of (1, 2, 3, 4, 5)." % atmtype 

658 raise Exception(errmsg) 

659 

660 # Inspect humidity (float). The range must be 0.0 gt. Humidity gt. 100.0 [%] 

661 humidity_float = float(humidity) 

662 if humidity_float != -1.0 and not (0.0 <= humidity_float and humidity_float <= 100.0): 

663 errmsg = "humidity (=%s) should be in range 0~100" % humidity 

664 raise Exception(errmsg) 

665 

666 # datacolumn check (by XML definition) 

667 datacolumn_upper = datacolumn.upper() 

668 if datacolumn_upper not in ['DATA', 'CORRECTED', 'FLOAT_DATA']: 

669 errmsg = "Specified column name (%s) Unacceptable." % datacolumn 

670 raise Exception(errmsg) 

671 

672 # tweak antenna selection string to include autocorr data 

673 antenna_autocorr = sdutil.get_antenna_selection_include_autocorr(infile, antenna) 

674 

675 # C++ re-implementation 

676 sdms.open(infile) 

677 sdms.set_selection(spw=outputspw, field=field, 

678 antenna=antenna_autocorr, 

679 timerange=timerange, scan=scan, 

680 polarization=correlation, intent=intent, 

681 observation=observation, feed=feed, 

682 taql=msselect, 

683 reindex=False) 

684 

685 config = get_configuration_for_atmcor( 

686 infile=infile, 

687 spw=spw, 

688 outputspw=outputspw, 

689 gainfactor=gainfactor, 

690 user_inputs=locals() 

691 ) 

692 

693 sdms.atmcor(config=config, datacolumn=datacolumn, outfile=outfile) 

694 

695 except Exception as err: 

696 casalog.post('%s' % err, priority='SEVERE') 

697 raise