Merge pull request #191 from 1054/20210320.vis.require.full.line.coverage

Check that cube frequency range is completely covered by all input Measurement Sets

Author: low-sky

phangsPipeline/casaVisRoutines.py

@@ -430,6 +430,7 @@ def find_spws_for_line(
     vsys_kms=None, vwidth_kms=None, vlow_kms=None, vhigh_kms=None,
     max_chanwidth_kms = None,
     require_data = False,
+    require_full_line_coverage = False, 
     exit_on_error = True, 
     as_list = False,
     ):
@@ -495,6 +496,8 @@ def find_spws_for_line(
     # SPWs one at a time.
 
     spw_list = []
+    spw_lowest_ghz = None
+    spw_highest_ghz = None
 
     logger.debug("... vm = au.ValueMapping(infile) ...")
     vm = au.ValueMapping(infile)
@@ -522,12 +525,23 @@ def find_spws_for_line(
                 continue
 
         spw_list.append(this_spw)
+        
+        if spw_lowest_ghz is None:
+            spw_lowest_ghz = spw_low_ghz
+        else:
+            spw_lowest_ghz = min(spw_lowest_ghz, spw_low_ghz)
+        
+        if spw_highest_ghz is None:
+            spw_highest_ghz = spw_high_ghz
+        else:
+            spw_highest_ghz = max(spw_highest_ghz, spw_high_ghz)
 
     # If we don't find the line in this data set, issue a warning and
     # return.
 
     if len(spw_list) == 0:
         logger.warning('No spectral windows contain the input line.')
+        spw_list = []
         spw_list_string = None # can't be '', that selects all
     else:
 
@@ -536,7 +550,16 @@ def find_spws_for_line(
 
         # make spw_list_string appropriate for use in selection
         spw_list_string = ','.join(np.array(spw_list).astype(str))
-     
+    
+    # check if the spws in this measurement set completely covers the given line frequency range
+    if not (spw_lowest_ghz <= line_low_ghz and spw_highest_ghz >= line_high_ghz):
+        logger.warning('The spectral windows in this measurement set "%s" (%.6f -- %.6f) does not cover the full "%s" line frequency range (%.6f -- %.6f).'%(\
+                        infile, spw_lowest_ghz, spw_highest_ghz, line, line_low_ghz, line_high_ghz))
+        if require_full_line_coverage:
+            spw_list = []
+            spw_list_string = None # can't be '', that selects all
+    
+    # return
     if as_list:
         return(spw_list)
     else:
@@ -689,6 +712,7 @@ def compute_common_chanwidth(
     vwidth_kms = None, 
     vlow_kms = None, 
     vhigh_kms = None, 
+    require_full_line_coverage = False, 
     ): 
     """
     Calculates the coarsest channel width among all spectral windows
@@ -720,7 +744,8 @@ def compute_common_chanwidth(
     coarsest_channel = None
     for this_infile in infile_list:
         # Find spws for line
-        spw_list_string = find_spws_for_line(this_infile, line, vsys_kms = vsys_kms, vwidth_kms = vwidth_kms)
+        spw_list_string = find_spws_for_line(this_infile, line, vsys_kms = vsys_kms, vwidth_kms = vwidth_kms, 
+            require_full_line_coverage = require_full_line_coverage)
 
         chan_widths_hz = au.getChanWidths(this_infile, spw_list_string)
 
@@ -751,7 +776,8 @@ def batch_extract_line(
     method = 'regrid_then_rebin',
     exact = False,
     freqtol = '',
-    clear_pointing = True,
+    clear_pointing = True, 
+    require_full_line_coverage = False, 
     overwrite = False,
     ):
     """
@@ -806,6 +832,7 @@ def batch_extract_line(
 
             this_scheme['outfile'] = this_outfile
             this_scheme['overwrite'] = overwrite
+            this_scheme['require_full_line_coverage'] = require_full_line_coverage
             split_file_list.append(this_outfile)
 
             # Execute line extraction
@@ -1042,6 +1069,7 @@ def extract_line(
     target_chan_kms = None,
     nchan = None,
     binfactor = None,
+    require_full_line_coverage = False, 
     overwrite = False, 
     ):
     """
@@ -1133,9 +1161,10 @@ def extract_line(
         spw_list = find_spws_for_line(
             infile = line, restfreq_ghz = restfreq_ghz,
             vsys_kms=vsys_kms, vwidth_kms=vwidth_kms, vlow_kms=vlow_kms, vhigh_kms=vhigh_kms,
-            require_data = True, exit_on_error = True, as_list = True,
+            require_data = True, require_full_line_coverage = require_full_line_coverage, 
+            exit_on_error = True, as_list = True,
             )
-        if spw_list is None:
+        if spw_list is None or len(spw_list) == 0:
             logging.error("No SPWs for selected line and velocity range.")
             return()
         spw = spw_list.join(',')
@@ -1161,7 +1190,9 @@ def extract_line(
             infile=infile, outfile=outfile, restfreq_ghz=restfreq_ghz, spw=spw,
             vstart_kms=vstart_kms, vwidth_kms=vwidth_kms,
             target_chan_kms=target_chan_kms, nchan=nchan, 
-            method='regrid')
+            method='regrid', 
+            require_full_line_coverage=require_full_line_coverage,
+            )
 
     if method == 'just_rebin' or method == 'regrid_then_rebin' or \
             method == 'rebin_then_regrid':
@@ -1170,12 +1201,16 @@ def extract_line(
 
         rebin_params, rebin_msg =  build_mstransform_call(
             infile=infile, outfile=outfile, restfreq_ghz=restfreq_ghz, spw=spw,
-            binfactor=binfactor, method='rebin')
+            binfactor=binfactor, method='rebin', 
+            require_full_line_coverage=require_full_line_coverage,
+            )
 
     if multiple_spws:
         combine_params, combine_msg =  build_mstransform_call(
             infile=infile, outfile=outfile, restfreq_ghz=restfreq_ghz, spw=spw,
-            method='combine')
+            method='combine', 
+            require_full_line_coverage=require_full_line_coverage,
+            )
 
     # ............................................
     # string the calls together in the desired order
@@ -1290,6 +1325,7 @@ def build_mstransform_call(
     target_chan_kms = None,
     nchan = None,
     binfactor = None,
+    require_full_line_coverage = False, 
     overwrite = False, 
     ):
     """    
@@ -1326,8 +1362,9 @@ def build_mstransform_call(
     if spw is None:
         if restfreq_ghz is not None:
             spw = find_spws_for_line(
-                infile=infile, restfreq_ghz = restfreq_ghz,
-                vlow_kms=vlow_kms, vhigh_kms=vhigh_kms)
+                infile=infile, restfreq_ghz=restfreq_ghz,
+                vlow_kms=vlow_kms, vhigh_kms=vhigh_kms,
+                require_full_line_coverage=require_full_line_coverage)
 
             # Exit if no SPWs contain the line.
             if spw is None:

phangsPipeline/handlerVis.py

@@ -101,15 +101,22 @@ def loop_stage_uvdata(
         statwt_cont = True,
         collapse_cont = True,
         timebin = None,
-        just_projects=None,
-        strict_config=True,
+        just_projects = None,
+        strict_config = True,
+        require_full_line_coverage = False,
         overwrite = False,
         ):
         """
         Loops over the full set of targets, products, and configurations
         to run the uv data processing. Toggle the parts of the loop
         using the do_XXX booleans. Other choices affect the algorithms
         used.
+        
+        The strict_config option sets whether to require that a target has data
+        from ALL arrays that make up the configuration (True) or not (False).
+        
+        The require_full_line_coverage option sets whether to require a measurement set 
+        to completely cover a given line's frequency range (True) or not (False). 
         """
 
         if make_directories:
@@ -153,7 +160,8 @@ def loop_stage_uvdata(
                             obsnum = this_obsnum,
                             product = this_product,
                             timebin = timebin,
-                            overwrite = overwrite,
+                            require_full_line_coverage = require_full_line_coverage, 
+                            overwrite = overwrite, 
                             )
 
                     # Run custom processing. Not currently used.
@@ -208,6 +216,7 @@ def loop_stage_uvdata(
                             extra_ext_in = "",
                             contsub = "prefer",
                             # could add algorithm flags here
+                            require_full_line_coverage = require_full_line_coverage, 
                             overwrite = overwrite,
                             )
 
@@ -278,6 +287,7 @@ def task_split(
             do_statwt = False,
             timebin = None,
             use_symlink = True,
+            require_full_line_coverage = False, 
             overwrite = False,
             ):
         """
@@ -366,7 +376,8 @@ def task_split(
                         spw = cvr.find_spws_for_line(
                             infile = infile, line = this_line,
                             max_chanwidth_kms = max_chanwidth_kms,
-                            vsys_kms = vsys, vwidth_kms = vwidth)
+                            vsys_kms = vsys, vwidth_kms = vwidth,
+                            require_full_line_coverage = require_full_line_coverage)
 
                     if spw is None or len(spw) == 0:
                         logger.warning("... No SPWs meet the selection criteria. Skipping.")
@@ -647,6 +658,7 @@ def task_extract_line(
             do_statwt = True,
             edge_for_statwt = None,
             method = "regrid_then_rebin",
+            require_full_line_coverage = False, 
             overwrite = False,
             ):
         """
@@ -803,6 +815,7 @@ def task_extract_line(
                 exact = exact,
                 overwrite = overwrite,
                 clear_pointing = False,
+                require_full_line_coverage = require_full_line_coverage, 
                 )
 
             if do_statwt: