Commits

Wataru Kawasaki authored 4b0eede2642
add unit tests for clipping
No tags

casatasks/tests/tasks/test_sdbaseline.py

Modified
6 6 import os
7 7 import shutil
8 8 import unittest
9 9
10 10 from casatasks import casalog
11 11 from casatasks import sdbaseline
12 12 from casatasks.private.sdutil import table_manager
13 13 from casatasks.private.task_sdbaseline import check_fftthresh
14 14 from casatasks.private.task_sdbaseline import is_empty
15 15 from casatasks.private.task_sdbaseline import parse_wavenumber_param
16 -# from casatestutils import selection_syntax
17 -from casatools import ctsys, table
16 +from casatools import ctsys
17 +from casatools import table
18 18
19 19
20 20 tb = table()
21 21 ctsys_resolve = ctsys.resolve
22 22
23 23
24 24 # default is necessary in CASA6
25 25 def default(atask):
26 26 pass
27 27
5241 5241
5242 5242 def test052(self):
5243 5243 self.params['spw'] = self.spw
5244 5244 self.run_apply_test()
5245 5245
5246 5246
5247 5247 class sdbaseline_clipping(sdbaseline_unittest_base):
5248 5248 """
5249 5249 Tests for iterative sigma clipping
5250 5250
5251 - test000 --- blfunc='poly'
5252 - test001 --- blfunc='cspline'
5253 - test002 --- blfunc='sinusoid'
5254 - test003 --- blfunc='variable'
5251 + test000 --- to confirm if clipping works regardless of blformat when blfunc='poly'
5252 + test001 --- to confirm if clipping works regardless of blformat blfunc='cspline'
5253 + test002 --- to confirm if clipping works regardless of blformat blfunc='sinusoid'
5254 + test003 --- to confirm if clipping works regardless of blformat blfunc='variable'
5255 +
5256 + test010 --- clipping runs multiple times (positive spikes only, threshold=3sigma)
5257 + test011 --- clipping runs multiple times (positive spikes only, threshold=10sigma)
5258 + test012 --- clipping runs multiple times (negative spikes only)
5259 + test013 --- clipping runs multiple times (both positive/negative spikes)
5260 +
5261 + test020 --- clipping does run but actually no data clipped (huge threshold)
5262 + test021 --- clipping does run but actually no data clipped (no spike)
5255 5263 """
5256 5264
5257 5265 datapath = ctsys_resolve('unittest/sdbaseline/')
5258 5266 infile = 'analytic_order3_withoffset.ms'
5259 5267 outroot = sdbaseline_unittest_base.taskname + '_clippingtest'
5260 5268 outfile = outroot + '.ms'
5269 + csvfile = infile + '_blparam.csv'
5261 5270 blparamfile = outroot + '.blparam'
5262 5271 params = {'infile': infile, 'outfile': outfile,
5263 5272 'overwrite': True,
5264 5273 'intent': 'OBSERVE_TARGET#ON_SOURCE',
5265 5274 'datacolumn': 'float_data',
5266 5275 'blmode': 'fit',
5267 5276 'order': 0, 'npiece': 1, 'addwn': 0,
5268 - 'clipthresh': 3.0,
5269 5277 'blparam': blparamfile}
5270 5278 outdata = {}
5279 + outmask = {}
5271 5280
5272 - def _setup_input_data(self):
5281 + def _setup_input_data(self, spikes):
5273 5282 with table_manager(self.infile, nomodify=False) as tb:
5274 5283 data = tb.getcell('FLOAT_DATA', 0)
5275 5284 for ipol in range(len(data)):
5285 + # base data by repeating 1 and -1 - this has mean=0, sigma=1
5276 5286 for ichan in range(len(data[0])):
5277 5287 data[ipol][ichan] = 1.0 if ichan % 2 == 0 else -1.0
5278 - data[ipol][4000] = 100000.0
5288 + # add spike data
5289 + for chan, value in spikes:
5290 + data[ipol][chan] = value
5279 5291 tb.putcell('FLOAT_DATA', 0, data)
5280 5292
5281 - def _set_params(self, blfunc, outbl, doclip):
5293 + def _set_params(self, blfunc, outbl, clipniter, thres):
5282 5294 self.params['blfunc'] = blfunc
5283 5295 self.params['blformat'] = 'csv' if outbl else ''
5284 - self.params['clipniter'] = 1 if doclip else 0
5296 + self.params['clipniter'] = clipniter
5297 + self.params['clipthresh'] = thres
5285 5298
5286 - def _get_data_name(self, outbl, doclip):
5299 + def _get_data_name(self, outbl, clipniter):
5287 5300 name_bl = 'bl' if outbl else 'nobl'
5288 - name_cl = 'clip' if doclip else 'noclip'
5301 + name_cl = str(clipniter)
5289 5302
5290 5303 return name_bl + '-' + name_cl
5291 5304
5292 - def _exec_sdbaseline(self, blfunc, outbl, doclip):
5293 - self._set_params(blfunc, outbl, doclip)
5305 + def _exec_sdbaseline(self, blfunc, outbl, clipniter, thres):
5306 + self._set_params(blfunc, outbl, clipniter, thres)
5294 5307
5295 5308 if blfunc == 'variable':
5296 - self._create_blparam_file(doclip)
5309 + self._create_blparam_file(clipniter, thres)
5297 5310
5298 5311 sdbaseline(**self.params)
5299 5312
5300 5313 with table_manager(self.outfile) as tb:
5301 - data_name = self._get_data_name(outbl, doclip)
5302 - # use data at row 0, pol 0 only
5303 - self.outdata[data_name] = tb.getcell('FLOAT_DATA', 0)[0]
5314 + data_name = self._get_data_name(outbl, clipniter)
5315 + self.outdata[data_name] = tb.getcell('FLOAT_DATA', 0)[0] # row(spw)=0, pol=0
5316 + if outbl:
5317 + with open(self.csvfile) as f:
5318 + for line in csv.reader(f):
5319 + if (line[2] == '0') and (line[3] == '0'): # row(spw)=0, pol=0
5320 + data_name = self._get_data_name(outbl, clipniter)
5321 + self.outmask[data_name] = line[5]
5322 +
5323 + remove_single_file_dir(self.csvfile)
5304 5324
5305 5325 remove_files_dirs(self.outroot)
5306 5326
5307 - def _result(self, outbl, doclip):
5308 - return self.outdata[self._get_data_name(outbl, doclip)]
5327 + def _result(self, outbl, clipniter):
5328 + return self.outdata[self._get_data_name(outbl, clipniter)]
5309 5329
5310 - def _create_blparam_file(self, doclip):
5311 - clipniter = str(1 if doclip else 0)
5330 + def _resmask(self, clipniter):
5331 + return self.outmask[self._get_data_name(True, clipniter)]
5312 5332
5333 + def _create_blparam_file(self, clipniter, thres):
5313 5334 with open(self.blparamfile, 'w') as f:
5314 - f.write('0,0,,' + clipniter + ',3.,false,,,,,poly,0,,[]')
5315 - f.write('0,1,,' + clipniter + ',3.,false,,,,,chebyshev,0,,[]')
5316 - f.write('1,0,,' + clipniter + ',3.,false,,,,,cspline,,1,[]')
5317 - f.write('1,1,,' + clipniter + ',3.,false,,,,,cspline,,1,[]')
5335 + f.write(f'0,0,,{clipniter},{thres},false,,,,,poly,0,,[]')
5336 + f.write(f'0,1,,{clipniter},{thres},false,,,,,chebyshev,0,,[]')
5337 + f.write(f'1,0,,{clipniter},{thres},false,,,,,cspline,,1,[]')
5338 + f.write(f'1,1,,{clipniter},{thres},false,,,,,cspline,,1,[]')
5339 +
5340 + def _run_test(self, blfunc='poly', spikes=[(4000, 100000.0)], thres=3.0, ifclipped=True):
5341 + self._setup_input_data(spikes=spikes)
5318 5342
5319 - def _run_test(self, blfunc):
5320 5343 bools = [False, True]
5321 - lst = [(outbl, doclip) for outbl in bools for doclip in bools]
5322 - for outbl, doclip in lst:
5323 - self._exec_sdbaseline(blfunc, outbl, doclip)
5344 + lst = [(outbl, clipniter) for outbl in bools for clipniter in [0, 1]]
5345 + for outbl, clipniter in lst:
5346 + self._exec_sdbaseline(blfunc, outbl, clipniter, thres)
5324 5347
5325 5348 # if clipping is turned on, output of sdbaseline must be identical
5326 5349 # regardless of whether blformat is empty or not
5327 - self.assertTrue(np.allclose(self._result(False, True), self._result(True, True)),
5350 + self.assertTrue(np.allclose(self._result(False, 1), self._result(True, 1)),
5328 5351 msg='unexpected result; result differs with different blformat.')
5329 5352 # with iterative clipping, output of sdbaseline must be different from that
5330 5353 # without clipping, regardless of whether blformat is empty or not
5331 - for blout in bools:
5332 - self.assertFalse(np.allclose(self._result(blout, True), self._result(blout, False)),
5333 - msg='unexpected result; clipping is not working.')
5354 + if ifclipped:
5355 + for blout in bools:
5356 + self.assertFalse(np.allclose(self._result(blout, 1), self._result(blout, 0)),
5357 + msg='unexpected result; clipping is not working.')
5358 + else:
5359 + for blout in bools:
5360 + self.assertTrue(np.allclose(self._result(blout, 1), self._result(blout, 0)),
5361 + msg='unexpected result; clipping is done.')
5362 +
5363 + def _run_test_multiple_clipping(self, spikes, thres=3.0):
5364 + # using a data with two spikes with different values.
5365 + # this test is to confirm if the larger spike is clipped in the first turn
5366 + # and if the smaller spike is clipped in the second turn
5367 + # and if no more data is clipped afterwards.
5368 +
5369 + self._setup_input_data(spikes=[(2000, 1000.0), (4000, 100000.0)])
5370 + maxiter = 5
5371 + for niter in range(maxiter):
5372 + self._exec_sdbaseline(blfunc='poly', outbl=True, clipniter=niter, thres=thres)
5373 +
5374 + answer_mask_before = '[[0;8191]]'
5375 + self.assertEqual(self._resmask(0), answer_mask_before)
5376 + answer_mask_after_1clip = '[[0;3999];[4001;8191]]' # channel 4000 is clipped
5377 + self.assertEqual(self._resmask(1), answer_mask_after_1clip)
5378 + answer_mask_after_2clip = '[[0;1999];[2001;3999];[4001;8191]]' # channel 2000 is clipped
5379 + for i in range(2, maxiter):
5380 + self.assertEqual(self._resmask(i), answer_mask_after_2clip)
5334 5381
5335 5382 def setUp(self):
5336 5383 remove_files_dirs(self.infile)
5337 5384 shutil.copytree(os.path.join(self.datapath, self.infile), self.infile)
5338 - self._setup_input_data()
5339 5385 default(sdbaseline)
5386 + self.outdata = {}
5387 + self.outmask = {}
5340 5388
5341 5389 def tearDown(self):
5342 5390 remove_files_dirs(self.infile)
5343 5391 remove_files_dirs(self.outroot)
5344 5392
5345 5393 def test000(self):
5346 - self._run_test('poly')
5394 + self._run_test(blfunc='poly')
5347 5395
5348 5396 def test001(self):
5349 - self._run_test('cspline')
5397 + self._run_test(blfunc='cspline')
5350 5398
5351 5399 def test002(self):
5352 - self._run_test('sinusoid')
5400 + self._run_test(blfunc='sinusoid')
5353 5401
5354 5402 def test003(self):
5355 - self._run_test('variable')
5403 + self._run_test(blfunc='variable')
5404 +
5405 + def test010(self):
5406 + self._run_test_multiple_clipping(spikes=[(2000, 1000.0), (4000, 100000.0)])
5407 +
5408 + def test011(self):
5409 + self._run_test_multiple_clipping(spikes=[(2000, 1000.0), (4000, 100000.0)], thres=10.0)
5410 +
5411 + def test012(self):
5412 + self._run_test_multiple_clipping(spikes=[(2000, -1000.0), (4000, -100000.0)])
5413 +
5414 + def test013(self):
5415 + self._run_test_multiple_clipping(spikes=[(2000, -1000.0), (4000, 100000.0)])
5416 +
5417 + def test020(self):
5418 + self._run_test(thres=100.0, ifclipped=False)
5419 +
5420 + def test021(self):
5421 + self._run_test(thres=3.0, spikes=[], ifclipped=False)
5356 5422
5357 5423
5358 5424 class sdbaseline_helperTest(sdbaseline_unittest_base):
5359 5425 """
5360 5426 Tests for helper functions
5361 5427
5362 5428 test000 --- tests for is_empty()
5363 5429 test010 --- tests for parse_wavenumber_param()
5364 5430 test020 --- tests for check_fftthresh()
5365 5431 """

Everything looks good. We'll let you know here if there's anything you should know about.

Add shortcut