casa / casa6 / 53c82146ccb

Commits

Takeshi Nakazato authored and Ville Suoranta committed 53c82146ccb Merge 01 Jul 2024

Pull request #828: CAS-14345

Merge in CASA/casa6 from CAS-14345 to master

* commit 'd1acbde20e5f9c8232703f5535ff6a2aa47b4169':
  Fix docstring issue
  using regex partial match (m) instead of full string match (p)
  Use TaQL to manipulate MS
  use onlyparse mode
  Fixed typo
  CAS-14345 test correction factor directly
  CAS-14345 add unit test for the fix
  CAS-14345 mitigate boundary effect of the convolution operation

casatasks/tests/tasks/test_task_sdatmcor.py

Modified

-27
+import os
-28
+import re
-29
+import shutil
-30
+import unittest
-31
-32
+import numpy as np
-33
-34
+from casatasks import applycal, casalog, gencal, sdatmcor
-35
+from casatasks.private.sdutil import (convert_antenna_spec_autocorr,
-36
+                                      get_antenna_selection_include_autocorr,
--
+                                      table_manager)
-+
+                                      table_manager, table_selector)
-38
+import casatasks.private.task_sdatmcor as sdatmcor_impl
-39
+from casatools import calibrater, ctsys
-40
+from casatools import ms as mstool
-41
+from casatools import quanta
-42
-43
+ctsys_resolve = ctsys.resolve
-44
-45
+# hold omp_get_num_threads value when this file is imported
-46
+OMP_NUM_THREADS_INITIAL = casalog.ompGetNumThreads()
-47
-644
+        casalog.post(
-645
+            f'OMP_NUM_THREAD_VALUES: initial: {OMP_NUM_THREADS_INITIAL} (returned by get_omp_num_threads), '
-646
+            f'at test start time: {omp_num_threads_org}, expected: {num_threads_expected}, '
-647
+            f'last set in logfile: {num_threads_log}')
-648
+        self.assertEqual(num_threads_expected, num_threads_log)
-649
-650
+        # check output MS
-651
+        self.check_result({19: True, 23: True})
-652
-653
-+
+def set_data_to_zero(infile: str, spw: int) -> np.ndarray:
-+
+    """Set ON_SOURCE CORRECTED_DATA to zero for given spw.
-+
-+
+    This is implemented by TaQL and is intended to do the
-+
+    same with the following code snippet.
-+
-+
+        ms.open(infile)
-+
+        ms.msselect(
-+
+            {'spw': str(spw), 'scanintent': 'OBSERVE_TARGET#ON_SOURCE'},
-+
+            onlyparse=True
-+
+        )
-+
+        msidx = ms.msselectedindices()
-+
+        ms.close()
-+
+        ddid = msidx['spwdd'][0]
-+
+        stateid = list(msidx['stateid'])
-+
+        taql = f'DATA_DESC_ID == {ddid} && STATE_ID IN {stateid}'
-+
-+
+        with table_selector(infile, taql=taql, nomodify=False) as tb:
-+
+            cdata = tb.getcol('CORRECTED_DATA')
-+
+            cdata[::] = 0
-+
+            tb.putcol(colname, cdata)
-+
+        return cdata
-+
-+
+    Args:
-+
+        infile: Input MS name
-+
+        spw: Spectral window Id
-+
-+
+    Returns:
-+
+        Data array manipulated by this function
-+
+    """
-+
+    with table_manager(infile) as tb:
-+
+        taql_string = f'''
-+
+        USING STYLE PYTHON
-+
+        UPDATE "{infile}" SET CORRECTED_DATA = 0
-+
+        WHERE
-+
+          DATA_DESC_ID IN
-+
+            [SELECT ROWID() FROM ::DATA_DESCRIPTION WHERE SPECTRAL_WINDOW_ID == {spw}]
-+
+          && STATE_ID IN
-+
+            [SELECT ROWID() FROM ::STATE WHERE OBS_MODE ~ m/^OBSERVE_TARGET#ON_SOURCE/]
-+
+        '''
-+
+        t = tb.taql(taql_string)
-+
+        cdata = t.getcol('CORRECTED_DATA')
-+
+        t.close()
-+
+    return cdata.real
-+
-+
-+
+class test_sdatmcor_smoothing(unittest.TestCase):
-+
+    datapath = 'measurementset/almasd'
-+
+    infile = 'X59ca_sel.ms'
-+
+    outfile = infile + '.atmcor'
-+
-+
+    def setUp(self):
-+
+        smart_remove(self.infile)
-+
+        smart_remove(self.outfile)
-+
+        datapath_to_ms = ctsys_resolve(os.path.join(self.datapath, self.infile))
-+
+        shutil.copytree(datapath_to_ms, self.infile)
-+
-+
+    def tearDown(self):
-+
+        smart_remove(self.infile)
-+
+        smart_remove(self.outfile)
-+
-+
+    def _get_data(self, ms_name, spw):
-+
+        with table_selector(ms_name, f'DATA_DESC_ID=={spw}') as tb:
-+
+            data = tb.getcol('DATA').real
-+
-+
+        return data
-+
-+
+    def test_mitigation(self):
-+
+        """Test if mitigation for boundary effect of convolution works."""
-+
+        # Set data zero to get correction factor
-+
+        # In sdatmcor, correction factor is subtracted from the data
-+
+        # so we can get correction factor if we set data all zero
-+
+        # (output_data = 0 - correction_factor).
-+
+        zero_data = set_data_to_zero(self.infile, spw=17)
-+
+        self.assertTrue(np.all(zero_data == 0))
-+
-+
+        # Apply correction with smoothing
-+
+        sdatmcor(
-+
+            infile=self.infile,
-+
+            outfile=self.outfile,
-+
+            spw='17,19',
-+
+            intent='OBSERVE_TARGET#ON_SOURCE',
-+
+            datacolumn='corrected',
-+
+            gainfactor={17: 41.49, 19: 41.48},
-+
+            dtem_dh=-5.6,
-+
+            h0=2.0,
-+
+            atmtype=1
-+
+        )
-+
-+
+        # Resulting data should be -correction_factor. Only check spectral
-+
+        # data for spw 17 since it is severely suffered from the boundary
-+
+        # effect. If mitigation didn't work, there will be steep increase
-+
+        # or decrease at edge channels. In that case, the slope should be
-+
+        # order of magnitude larger.
-+
+        correction_factor_spw17 = self._get_data(self.outfile, spw=17)
-+
+        average_factor_per_pol = correction_factor_spw17.mean(axis=2)
-+
+        # average derivative excluding edge channels
-+
+        delta = average_factor_per_pol[:, 1:] - average_factor_per_pol[:, :-1]
-+
+        average_delta = delta[:, 10:-10].mean()
-+
+        threshold = abs(average_delta) * 10
-+
+        for ipol in range(correction_factor_spw17.shape[0]):
-+
+            print(f'Examining pol {ipol}')
-+
+            edge_delta = np.abs(delta[ipol, [0, -1]])
-+
+            print('edge_delta', edge_delta)
-+
+            print(f'threshold = {threshold}')
-+
+            self.assertTrue(
-+
+                np.all(edge_delta < threshold),
-+
+                msg=f'Mitigation did not work for pol {ipol}'
-+
+            )
-+
-+
-765
+class ATMParamTest(unittest.TestCase):
-766
+    def _param_test_template(self, valid_test_cases,
-767
+                             invalid_user_input, user_default, task_default, unit=''):
-768
+        # internal error
-769
+        wrong_task_default = 'NG'
-770
+        with self.assertRaises(RuntimeError):
-771
+            param, is_customized = sdatmcor_impl.parse_atm_params(
-772
+                '',
-773
+                user_default,
-774
+                wrong_task_default

casa6

Commits

Jira Issues

Add shortcut