More restrictive conditions to detect weird oscillations in GSFC reco…

…rds (in order to be deleted); th_runDetect modified
bcobo · Feb 11, 2021 · 1c5f8aa · 1c5f8aa
1 parent cb5f46e
commit 1c5f8aa
Show file tree

Hide file tree

Showing 5 changed files with 46 additions and 13 deletions.
diff --git a/doc/SIRENAfunctions.rst b/doc/SIRENAfunctions.rst
@@ -5616,7 +5616,7 @@ Search functions by name at :ref:`genindex`.
         Pulseheight of the input pulse    
 
 
-.. cpp:function:: void runDetect(TesRecord* record, int trig_reclength, int lastRecord, PulsesCollection *pulsesAll, ReconstructInitSIRENA** reconstruct_init, PulsesCollection** pulsesInRecord)
+.. cpp:function:: void runDetect(TesRecord* record, int trig_reclength, int lastRecord, int nrecord, PulsesCollection *pulsesAll, ReconstructInitSIRENA** reconstruct_init, PulsesCollection** pulsesInRecord)
      
     Located in file: *tasksSIRENA.cpp*
 
@@ -5643,10 +5643,12 @@ Search functions by name at :ref:`genindex`.
         3) (Filter and) differentiate the *models* of the library (only for the first record in :option:`opmode` = 1). Run  (:cpp:func:`filderLibrary`)
 
         4) Store the input record in *invector* (:cpp:func:`loadRecord`)
+
+        5) Detect weird oscillations in some GSFC records
 
-        5) Convert *I* into *R* if :option:`EnergyMethod` = **I2R** or **I2RFITTED** (:cpp:func:`convertI2R`)
+        6) Convert *I* into *R* if :option:`EnergyMethod` = **I2R** or **I2RFITTED** (:cpp:func:`convertI2R`)
 
-        6) Process each record (:cpp:func:`proceRecord`): 
+        7) Process each record (:cpp:func:`proceRecord`): 
 
                 - (Low-pass filter and) differentiate                          
                 - Find pulses                                                           
@@ -5656,12 +5658,12 @@ Search functions by name at :ref:`genindex`.
 
         **From this point forward, I2R and I2RFITTED are completely equivalent to OPTFILT**     
 
-        7) If last record in :option:`opmode` = 0 run:                                         
+        8) If last record in :option:`opmode` = 0 run:                                         
 
                 * :cpp:func:`calculateTemplate` (and :cpp:func:`weightMatrix`)
                 * :cpp:func:`writeLibrary`
 
-        8) If last record and PCA:
+        9) If last record and PCA:
 
                 - In order to not have restrictions when providing (\*reconstruct_init)->energyPCAx
                 - Covariance data
@@ -5674,7 +5676,7 @@ Search functions by name at :ref:`genindex`.
                 - Conversion factor from arbitrary unit to eV
                 - Energy calculation
 
-        9) Close intermediate output FITS file if it is necessary   
+        10) Close intermediate output FITS file if it is necessary   
 
     **Members/Variables**
 
@@ -5714,6 +5716,10 @@ Search functions by name at :ref:`genindex`.
 
         Integer to verify whether *record* is the last one (=1) to be read (and thus if library file will be created)
 
+    .. cpp:member:: int nrecord
+
+        Current record index (to know the particular record where there is a weird oscillation)
+
     .. cpp:member:: PulsesCollection* pulsesAll
 
         Member of *PulsesCollection* structure to successively store all the pulses used to create the library. Re-populated after each processed record
@@ -6400,7 +6406,7 @@ Search functions by name at :ref:`genindex`.
         XML input FITS file with instrument definition
 
 
-.. cpp:function:: int th_runDetect (TesRecord* record, int trig_reclength, int lastRecord, PulsesCollection *pulsesAll, ReconstructInitSIRENA** reconstruct_init, PulsesCollection** pulsesInRecord)
+.. cpp:function:: int th_runDetect (TesRecord* record, int trig_reclength, int lastRecord, int nrecord, PulsesCollection *pulsesAll, ReconstructInitSIRENA** reconstruct_init, PulsesCollection** pulsesInRecord)
      
     Located in file: *tasksSIRENA.cpp*
 

diff --git a/libsixt/scheduler.cpp b/libsixt/scheduler.cpp
@@ -36,6 +36,7 @@ void detection_worker()
       //log_trace("Extracting detection data from queue...");
       th_runDetect(data->rec, data->trig_reclength,
                 data->last_record,
+                data->n_record,
                 data->all_pulses,
                 &(data->rec_init),
                 &(data->record_pulses));

diff --git a/libsixt/tasksSIRENA.cpp b/libsixt/tasksSIRENA.cpp
@@ -106,6 +106,7 @@
  * - Create intermediate output FITS file if required ('createDetectFile')
  * - (Filter and) differentiate the 'models' of the library (only for the first record in PRODUCTION 'opmode=1') ('filderLibrary')
  * - Store the input record in 'invector' ('loadRecord')
+ * - Detect weird oscillations in some GSFC records
  * - Convert I into R if 'EnergyMethod' = I2R or I2RFITTED ('convertI2R')
  * - Process each record ('proceRecord')
  * 	- (Low-pass filter and) differentiate
@@ -135,7 +136,7 @@
  * - record: Member of TesRecord' structure that contains the input record
  * - trig_reclength: Record size (just in case threading and input files with different 'ADC' lengths but the same record size indeed)
  * - lastRecord: Integer to verify whether record is the last one (=1) to be read (and thus if library file will be created)
- * - nrecord: Current record index
+ * - nrecord: Current record index (to know the particular record where there is a weird oscillation)
  * - pulsesAll: Member of 'PulsesCollection' structure to successively store all the pulses used to create the library. Re-populated after each processed record
  * - reconstruct_init: Member of 'ReconstructInitSIRENA' structure to initialize the reconstruction parameters (pointer and values)
  * - pulsesInRecord: Member of 'PulsesCollection' structure to store all the pulses found in the input record
@@ -239,9 +240,16 @@ void runDetect(TesRecord* record, int trig_reclength, int lastRecord, int nrecor
         message = "Cannot run findMeanSigma in runDetect";
         EP_EXIT_ERROR(message,EPFAIL);
     }
-    // sgTES~0.1% of meanTEST if the record is noise only => sgTES~2% of meanTEST will be weird oscillations but not noise
+    /*cout<<"meanTEST: "<<meanTEST<<endl;
+    cout<<"sgTEST: "<<sgTEST<<endl;
+    cout<<"max: "<<gsl_vector_max(invector)<<endl;
+    cout<<"min: "<<gsl_vector_min(invector)<<endl;
+    cout<<"gsl_vector_max(invector)-meanTEST: "<<gsl_vector_max(invector)-meanTEST<<endl;
+    cout<<"meanTEST-gsl_vector_min(invector): "<<meanTEST-gsl_vector_min(invector)<<endl;
+    cout<<"sgTEST*100/meanTEST: "<<sgTEST*100/meanTEST<<endl;*/
+    // sgTES~0.1% of meanTEST if the record is noise only => sgTES>1% is a weird oscillation
     // Noise records are not important if no pulses are detected in them
-    if ((gsl_vector_max(invector)-meanTEST) < (meanTEST-gsl_vector_min(invector)) && (gsl_vector_max(invector)-meanTEST>2*sgTEST) && (meanTEST-gsl_vector_min(invector)>2*sgTEST) && (sgTEST*100/meanTEST>2))
+    if ((gsl_vector_max(invector)-meanTEST) < (meanTEST-gsl_vector_min(invector)) && ((gsl_vector_max(invector)-meanTEST>sgTEST) || (meanTEST-gsl_vector_min(invector)>sgTEST)) && (meanTEST-gsl_vector_min(invector)>2*sgTEST) && (sgTEST*100/meanTEST>1))
     {
         char str_nrecord[125];      snprintf(str_nrecord,125,"%d",nrecord);
         message = "Weird oscillations in record " + string(str_nrecord);
@@ -780,7 +788,7 @@ void runDetect(TesRecord* record, int trig_reclength, int lastRecord, int nrecor
 std::mutex library_mut;
 std::mutex fits_file_mut;
 
-void th_runDetect(TesRecord* record, int trig_reclength, int lastRecord, PulsesCollection *pulsesAll, ReconstructInitSIRENA** reconstruct_init, PulsesCollection** pulsesInRecord)
+void th_runDetect(TesRecord* record, int trig_reclength, int lastRecord, int nrecord, PulsesCollection *pulsesAll, ReconstructInitSIRENA** reconstruct_init, PulsesCollection** pulsesInRecord)
 {
     //log_trace("th_runDetect: START");
     scheduler* sc = scheduler::get();
@@ -854,6 +862,23 @@ void th_runDetect(TesRecord* record, int trig_reclength, int lastRecord, PulsesC
     gsl_vector *invectorOriginal = gsl_vector_alloc(invector->size);
     gsl_vector_memcpy(invectorOriginal,invector);
 
+     // To detect weird oscillations in some GSFC records
+    double meanTEST=0;
+    double sgTEST=0;
+    if (findMeanSigma (invector, &meanTEST, &sgTEST))
+    {
+        message = "Cannot run findMeanSigma in runDetect";
+        EP_EXIT_ERROR(message,EPFAIL);
+    }
+    // sgTES~0.1% of meanTEST if the record is noise only => sgTES>1% is a weird oscillation
+    // Noise records are not important if no pulses are detected in them
+    if ((gsl_vector_max(invector)-meanTEST) < (meanTEST-gsl_vector_min(invector)) && ((gsl_vector_max(invector)-meanTEST>sgTEST) || (meanTEST-gsl_vector_min(invector)>sgTEST)) && (meanTEST-gsl_vector_min(invector)>2*sgTEST) && (sgTEST*100/meanTEST>1))
+    {
+        char str_nrecord[125];      snprintf(str_nrecord,125,"%d",nrecord);
+        message = "Weird oscillations in record " + string(str_nrecord);
+        EP_EXIT_ERROR(message,EPFAIL);
+    }
+
     if ((strcmp((*reconstruct_init)->EnergyMethod,"I2R") == 0) || (strcmp((*reconstruct_init)->EnergyMethod,"I2RFITTED") == 0))
     {
         if (convertI2R((*reconstruct_init)->EnergyMethod,(*reconstruct_init)->i2rdata->I0_START,(*reconstruct_init)->i2rdata->IMIN,(*reconstruct_init)->i2rdata->IMAX,(*reconstruct_init)->i2rdata->ADU_CNV, (*reconstruct_init)->i2rdata->ADU_BIAS,(*reconstruct_init)->i2rdata->I_BIAS,(*reconstruct_init)->i2rdata->Ifit,1/record->delta_t,&invector))

diff --git a/libsixt/tasksSIRENA.h b/libsixt/tasksSIRENA.h
@@ -54,6 +54,7 @@ void runDetect(TesRecord* record,
 
 void th_runDetect(TesRecord* record, int trig_reclength,
                   int lastRecord, 
+                  int nrecord,
                   PulsesCollection *pulsesAll, 
                   ReconstructInitSIRENA** reconstruct_init, 
                   PulsesCollection** pulsesInRecord);

diff --git a/libsixt/versionSIRENA.h b/libsixt/versionSIRENA.h
@@ -22,11 +22,11 @@
 //   CANTABRIA (CSIC-UC) with funding from the Spanish Ministry of Science and 
 //   Innovation (MICINN)
 //
-//   DATE: 2021/02/10, 17:08:31
+//   DATE: 2021/02/11, 10:38:37
 
 #ifndef SIRENA_VERSION_H
 #define SIRENA_VERSION_H
 
-#define SIRENA_VERSION "3.8.4"
+#define SIRENA_VERSION "3.8.5"
 
 #endif