LCOV - code coverage report
Current view: top level - flagging/Flagging - FlagAgentExtension.cc (source / functions) Hit Total Coverage
Test: ctest_coverage.info Lines: 83 123 67.5 %
Date: 2023-11-06 10:06:49 Functions: 5 5 100.0 %

          Line data    Source code
       1             : //# FlagAgenExtension.cc: This file contains the implementation of the FlagAgenExtension.cc class.
       2             : //#
       3             : //#  CASA - Common Astronomy Software Applications (http://casa.nrao.edu/)
       4             : //#  Copyright (C) Associated Universities, Inc. Washington DC, USA 2011, All rights reserved.
       5             : //#  Copyright (C) European Southern Observatory, 2011, All rights reserved.
       6             : //#
       7             : //#  This library is free software; you can redistribute it and/or
       8             : //#  modify it under the terms of the GNU Lesser General Public
       9             : //#  License as published by the Free software Foundation; either
      10             : //#  version 2.1 of the License, or (at your option) any later version.
      11             : //#
      12             : //#  This library is distributed in the hope that it will be useful,
      13             : //#  but WITHOUT ANY WARRANTY, without even the implied warranty of
      14             : //#  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
      15             : //#  Lesser General Public License for more details.
      16             : //#
      17             : //#  You should have received a copy of the GNU Lesser General Public
      18             : //#  License along with this library; if not, write to the Free Software
      19             : //#  Foundation, Inc., 59 Temple Place, Suite 330, Boston,
      20             : //#  MA 02111-1307  USA
      21             : //# $Id: $
      22             : 
      23             : #include <flagging/Flagging/FlagAgentExtension.h>
      24             : 
      25             : using namespace casacore;
      26             : namespace casa { //# NAMESPACE CASA - BEGIN
      27             : 
      28          29 : FlagAgentExtension::FlagAgentExtension(FlagDataHandler *dh, Record config, Bool writePrivateFlagCube):
      29          29 :                 FlagAgentBase(dh,config,ANTENNA_PAIRS_FLAGS,writePrivateFlagCube)
      30             : {
      31          29 :         setAgentParameters(config);
      32             : 
      33             :         // Request loading polarization map to FlagDataHandler
      34          29 :         flagDataHandler_p->setMapPolarizations(true);
      35          29 : }
      36             : 
      37          58 : FlagAgentExtension::~FlagAgentExtension()
      38             : {
      39             :         // Compiler automagically calls FlagAgentBase::~FlagAgentBase()
      40          58 : }
      41             : 
      42          29 : void FlagAgentExtension::setAgentParameters(Record config)
      43             : {
      44          29 :         logger_p->origin(LogOrigin(agentName_p,__FUNCTION__,WHERE));
      45             :         int exists;
      46             : 
      47          29 :         exists = config.fieldNumber ("extendpols");
      48          29 :         if (exists >= 0)
      49             :         {
      50          29 :                 if( config.type(exists) != TpBool )
      51             :                 {
      52           0 :                          throw( AipsError ( "Parameter 'extendpols' must be of type 'bool'" ) );
      53             :                 }
      54             :                 
      55          29 :                 extendpols_p = config.asBool("extendpols");
      56             :         }
      57             :         else
      58             :         {
      59           0 :                 extendpols_p = true;
      60             :         }
      61             : 
      62          29 :         *logger_p << logLevel_p << " extendpols is " << extendpols_p << LogIO::POST;
      63             : 
      64             : 
      65          29 :         exists = config.fieldNumber ("growtime");
      66          29 :         if (exists >= 0)
      67             :         {
      68          29 :                 if( config.type(exists) != TpDouble && config.type(exists) != TpFloat  && config.type(exists) != TpInt)
      69             :                 {
      70           0 :                          throw( AipsError ( "Parameter 'growtime' must be of type 'double'" ) );
      71             :                 }
      72             :                 
      73          29 :                 growtime_p = config.asDouble("growtime");
      74             : 
      75          29 :                 if( growtime_p < 0.0 || growtime_p > 100.0 )
      76             :                 {
      77           0 :                          throw( AipsError ( "Unsupported value for growtime:" + String::toString(growtime_p) + ". Allowed range is 0.0 through 100.0" ) );
      78             :                 }
      79             :         }
      80             :         else
      81             :         {
      82           0 :                 growtime_p = 50.0;
      83             :         }
      84             : 
      85          29 :         *logger_p << logLevel_p << " growtime is " << growtime_p << LogIO::POST;
      86             : 
      87             : 
      88          29 :         exists = config.fieldNumber ("growfreq");
      89          29 :         if (exists >= 0)
      90             :         {
      91          29 :                 if( config.type(exists) != TpDouble && config.type(exists) != TpFloat  && config.type(exists) != TpInt )
      92             :                 {
      93           0 :                          throw( AipsError ( "Parameter 'growfreq' must be of type 'double'" ) );
      94             :                 }
      95             :                 
      96          29 :                 growfreq_p = config.asDouble("growfreq");
      97             : 
      98          29 :                 if( growfreq_p < 0.0 || growfreq_p > 100.0 )
      99             :                 {
     100           0 :                          throw( AipsError ( "Unsupported value for growfreq:" + String::toString(growfreq_p) + ". Allowed range is 0.0 through 100.0" ) );
     101             :                 }
     102             :         }
     103             :         else
     104             :         {
     105           0 :                 growfreq_p = 50.0;
     106             :         }
     107             : 
     108          29 :         *logger_p << logLevel_p << " growfreq is " << growfreq_p << LogIO::POST;
     109             : 
     110             : 
     111          29 :         exists = config.fieldNumber ("growaround");
     112          29 :         if (exists >= 0)
     113             :         {
     114          27 :                 if( config.type(exists) != TpBool )
     115             :                 {
     116           0 :                          throw( AipsError ( "Parameter 'growaround' must be of type 'bool'" ) );
     117             :                 }
     118             :                 
     119          27 :                 growaround_p = config.asBool("growaround");
     120             :         }
     121             :         else
     122             :         {
     123           2 :                 growaround_p = false;
     124             :         }
     125             : 
     126          29 :         *logger_p << logLevel_p << " growaround is " << growaround_p << LogIO::POST;
     127             : 
     128             : 
     129          29 :         exists = config.fieldNumber ("flagneartime");
     130          29 :         if (exists >= 0)
     131             :         {
     132          27 :                 if( config.type(exists) != TpBool )
     133             :                 {
     134           0 :                          throw( AipsError ( "Parameter 'flagneartime' must be of type 'bool'" ) );
     135             :                 }
     136             :                 
     137          27 :                 flagneartime_p = config.asBool("flagneartime");
     138             :         }
     139             :         else
     140             :         {
     141           2 :                 flagneartime_p = false;
     142             :         }
     143             : 
     144          29 :         *logger_p << logLevel_p << " flagneartime is " << flagneartime_p << LogIO::POST;
     145             : 
     146             : 
     147          29 :         exists = config.fieldNumber ("flagnearfreq");
     148          29 :         if (exists >= 0)
     149             :         {
     150          27 :                 if( config.type(exists) != TpBool )
     151             :                 {
     152           0 :                          throw( AipsError ( "Parameter 'flagnearfreq' must be of type 'bool'" ) );
     153             :                 }
     154             :                 
     155          27 :                 flagnearfreq_p = config.asBool("flagnearfreq");
     156             :         }
     157             :         else
     158             :         {
     159           2 :                 flagnearfreq_p = false;
     160             :         }
     161             : 
     162          29 :         *logger_p << logLevel_p << " flagnearfreq is " << flagnearfreq_p << LogIO::POST;
     163             : 
     164             : 
     165          29 :         return;
     166             : }
     167             : 
     168             : bool
     169        4812 : FlagAgentExtension::computeAntennaPairFlags(const vi::VisBuffer2 &/*visBuffer*/,
     170             :                                 FlagMapper &flags,Int /*antenna1*/,Int /*antenna2*/,vector<uInt> &/*rows*/)
     171             : {
     172             :         // Set logger origin
     173        4812 :         logger_p->origin(LogOrigin(agentName_p,__FUNCTION__,WHERE));
     174             : 
     175             :         // Get flag cube size
     176        4812 :         IPosition flagCubeShape = flags.shape();
     177             :         Int nPols,nChannels,nTimesteps;
     178        4812 :         flags.shape(nPols, nChannels, nTimesteps);
     179             : 
     180             :         // STEP1: Grow flags around, across time and across frequency
     181             :         bool flag;
     182             :         uInt extendAround;
     183        4812 :         if (flagneartime_p or flagnearfreq_p or growaround_p)
     184             :         {
     185             :                 // Create the growing flag cube setting all flags initially to false
     186           0 :                 Cube<Bool> *growFlags = new Cube<Bool>(nPols, nChannels+2, nTimesteps+2,false);
     187             : 
     188             :                 // Iterate over flag cube and grow flags
     189           0 :                 for (Int timestep_i=0;timestep_i<nTimesteps;timestep_i++)
     190             :                 {
     191           0 :                         for (Int chan_j=0;chan_j<nChannels;chan_j++)
     192             :                         {
     193           0 :                                 for (Int pol_k=0;pol_k<nPols;pol_k++)
     194             :                                 {
     195           0 :                                         flag = flags.getModifiedFlags(pol_k,chan_j,timestep_i);
     196             : 
     197             :                                         // If the point is flagged we may extend in the time,frequency or polarization direction
     198           0 :                                         if (flag)
     199             :                                         {
     200           0 :                                                 if (flagneartime_p)
     201             :                                                 {
     202             :                                                         // Extend flags to the next timestep
     203           0 :                                                         growFlags->operator()(pol_k,chan_j+1,timestep_i + 2) = true;
     204             : 
     205             :                                                         // Extend flags to the previous timestep
     206           0 :                                                         growFlags->operator()(pol_k,chan_j+1,timestep_i) = true;
     207             :                                                 }
     208             : 
     209           0 :                                                 if (flagnearfreq_p)
     210             :                                                 {
     211             :                                                         // Extend flags to the next frequency
     212           0 :                                                         growFlags->operator()(pol_k,chan_j+2,timestep_i+1) = true;
     213             : 
     214             :                                                         // Extend flags to the previous frequency
     215           0 :                                                         growFlags->operator()(pol_k,chan_j,timestep_i+1) = true;
     216             :                                                 }
     217             :                                         }
     218             :                                         // Only if the point is not flagged, and grow around is activated, we check the nearest neighbours
     219           0 :                                         else if (growaround_p)
     220             :                                         {
     221           0 :                                                 extendAround = 0;
     222           0 :                                                 for (Int timestep_inner=timestep_i-1;timestep_inner<=timestep_i+1;timestep_inner++)
     223             :                                                 {
     224             :                                                         // Calculate number of flagged neighbours
     225           0 :                                                         for (Int chan_inner=chan_j-1;chan_inner<=chan_j+1;chan_inner++)
     226             :                                                         {
     227           0 :                                                                 if (    (timestep_inner>=0) and
     228           0 :                                                                                 (timestep_inner<nTimesteps) and
     229           0 :                                                                                 (chan_inner>=0) and
     230           0 :                                                                                 (chan_inner<nChannels))
     231             :                                                                 {
     232           0 :                                                                         extendAround += flags.getModifiedFlags(pol_k,chan_inner,timestep_inner);
     233             :                                                                 }
     234             :                                                         } // End grow around channel iteration
     235             :                                                 } // End grow around time iteration
     236             : 
     237             :                                                 // Extend flags if we have more than 4 flagged  neighbours
     238           0 :                                                 if (extendAround>4) growFlags->operator()(pol_k,chan_j+1,timestep_i+1) = true;
     239             : 
     240             :                                         } // End grow around
     241             : 
     242             :                                 } // For polarization
     243             :                         } // For channel
     244             :                 } // For timestep
     245             : 
     246             :                 // Apply grow flags to the common modified flag cube
     247           0 :                 for (uInt timestep_i=0;timestep_i<static_cast<uInt>(nTimesteps);timestep_i++)
     248             :                 {
     249           0 :                         for (uInt chan_j=0;chan_j<static_cast<uInt>(nChannels);chan_j++)
     250             :                         {
     251           0 :                                 for (uInt pol_k=0;pol_k<static_cast<uInt>(nPols);pol_k++)
     252             :                                 {
     253           0 :                                         if (growFlags->operator()(pol_k,chan_j+1,timestep_i+1))
     254             :                                         {
     255           0 :                                                 flags.setModifiedFlags(pol_k,chan_j,timestep_i);
     256           0 :                                                 visBufferFlags_p += 1;
     257             :                                         }
     258             :                                 }
     259             :                         }
     260             :                 }
     261             : 
     262             :                 // Delete growing flag cube
     263           0 :                 delete growFlags;
     264             : 
     265             :         } // End grow flags step
     266             : 
     267             : 
     268             :         // STEP2: Extend flags across time and frequency
     269        4812 :         if ( (growtime_p > 0) or (growfreq_p > 0) )
     270             :         {
     271             :                 // Per correlation accumulate over time and frequency
     272        4488 :                 Matrix<Double> acummTime(nPols,nTimesteps,0);
     273        4488 :                 Matrix<Double> acummFreq(nPols,nChannels,0);
     274      117739 :                 for (Int timestep_i=0;timestep_i<nTimesteps;timestep_i++)
     275             :                 {
     276     6543245 :                         for (Int chan_j=0;chan_j<nChannels;chan_j++)
     277             :                         {
     278    31933606 :                                 for (Int pol_k=0;pol_k<nPols;pol_k++)
     279             :                                 {
     280    25505856 :                                         flag = flags.getModifiedFlags(pol_k,chan_j,timestep_i);
     281    25505856 :                                         acummTime(pol_k,timestep_i) += flag;
     282    25505856 :                                         acummFreq(pol_k,chan_j) += flag;
     283             :                                 }
     284             :                         }
     285             :                 }
     286             : 
     287             :                 // Extend across frequency
     288        2244 :                 if (growfreq_p > 0)
     289             :                 {
     290      100363 :                         for (Int timestep_i=0;timestep_i<nTimesteps;timestep_i++)
     291             :                         {
     292      488490 :                                 for (Int pol_k=0;pol_k<nPols;pol_k++)
     293             :                                 {
     294      390179 :                                         if ((100.0*acummTime(pol_k,timestep_i)/nChannels) >= growfreq_p)
     295             :                                         {
     296             :                                                 // Flag all the channels for this polarization and timestep
     297      498067 :                                                 for (Int chan_j=0;chan_j<nChannels;chan_j++)
     298             :                                                 {
     299      489928 :                                                         flags.setModifiedFlags(pol_k,chan_j,timestep_i);
     300      489928 :                                                         visBufferFlags_p += 1;
     301             :                                                 }
     302             :                                         }
     303             :                                 }
     304             :                         }
     305             :                 }
     306             : 
     307             :                 // Extend across time
     308        2244 :                 if (growtime_p > 0)
     309             :                 {
     310      117105 :                         for (Int chan_j=0;chan_j<nChannels;chan_j++)
     311             :                         {
     312      465563 :                                 for (Int pol_k=0;pol_k<nPols;pol_k++)
     313             :                                 {
     314      350510 :                                         if ((100.0*acummFreq(pol_k,chan_j)/nTimesteps) >= growtime_p)
     315             :                                         {
     316             :                                                 // Flag all the timesteps for this polarization and channel
     317      672978 :                                                 for (Int timestep_i=0;timestep_i<nTimesteps;timestep_i++)
     318             :                                                 {
     319      573032 :                                                         flags.setModifiedFlags(pol_k,chan_j,timestep_i);
     320      573032 :                                                         visBufferFlags_p += 1;
     321             :                                                 }
     322             :                                         }
     323             :                                 }
     324             :                         }
     325             :                 }
     326             :         }
     327             : 
     328             :         // STEP3: Extend across polarizations
     329        4812 :         if (extendpols_p)
     330             :         {
     331       96313 :                 for (Int timestep_i=0;timestep_i<nTimesteps;timestep_i++)
     332             :                 {
     333     5073107 :                         for (Int chan_j=0;chan_j<nChannels;chan_j++)
     334             :                         {
     335    22377004 :                                 for (Int pol_k=0;pol_k<nPols;pol_k++)
     336             :                                 {
     337    17696904 :                                         if (flags.getModifiedFlags(pol_k,chan_j,timestep_i))
     338             :                                         {
     339     1250966 :                                                 for (Int ineer_pol_k=0;ineer_pol_k<nPols;ineer_pol_k++)
     340             :                                                 {
     341      949844 :                                                         flags.setModifiedFlags(ineer_pol_k,chan_j,timestep_i);
     342             :                                                 }
     343      301122 :                                                 visBufferFlags_p += nPols;
     344      301122 :                                                 break;
     345             :                                         }
     346             :                                 }
     347             :                         }
     348             :                 }
     349             :         }
     350             : 
     351        9624 :         return false;
     352             : }
     353             : 
     354             : } //# NAMESPACE CASA - END
     355             : 
     356             : 

Generated by: LCOV version 1.16