Filter Events (#239)

* Add FilterEvents tools
* Update README.md of FilterEvents toolchain
* removed FilterName config setting from FilterEvents tool
* add FilteredEvent tag to FilterEvents output file
* removed unnecessary Set(EventNumber)

UserTools/EventSelector/EventSelector.cpp

@@ -46,6 +46,7 @@ bool EventSelector::Initialise(std::string configfile, DataModel &data){
   m_variables.Get("SaveStatusToStore", fSaveStatusToStore);
   m_variables.Get("IsMC",fIsMC);
   m_variables.Get("RecoPDG",fRecoPDG);
+  m_variables.Get("CutConfiguration",fCutConfigurationName);
 
   if (!fIsMC){fMCFVCut = false; fMCPMTVolCut = false; fMCMRDCut = false; fMCPiKCut = false; fMCIsMuonCut = false; fMCIsElectronCut = false; fMCIsSingleRingCut = false; fMCIsMultiRingCut = false; fMCProjectedMRDHit = false; fMCEnergyCut = false; fPromptTrigOnly = false;}
 
@@ -69,6 +70,9 @@ bool EventSelector::Initialise(std::string configfile, DataModel &data){
 
   m_data->CStore.Get("ChannelNumToTankPMTSPEChargeMap",ChannelNumToTankPMTSPEChargeMap);
 
+  m_data->CStore.Set("CutConfiguration",fCutConfigurationName);
+
+
   return true;
 }
 

UserTools/EventSelector/EventSelector.h

@@ -237,6 +237,7 @@ class EventSelector: public Tool {
   bool fIsMC; 
   int fTriggerWord;
   int fRecoPDG;
+  std::string fCutConfigurationName;  
 
   bool get_mrd = false;
   double pmt_time = 0; 

UserTools/Factory/Factory.cpp

@@ -156,6 +156,7 @@ if (tool=="LAPPDDataDecoder") ret=new LAPPDDataDecoder;
 if (tool=="PythonScript") ret=new PythonScript;
 if (tool=="ReweightEventsGenie") ret=new ReweightEventsGenie;
 if (tool=="FilterLAPPDEvents") ret=new FilterLAPPDEvents;
+if (tool=="FilterEvents") ret=new FilterEvents;
 if (tool=="Stage1DataBuilder") ret=new Stage1DataBuilder;
 return ret;
 }

UserTools/FilterEvents/FilterEvents.cpp

@@ -0,0 +1,140 @@
+#include "FilterEvents.h"
+
+FilterEvents::FilterEvents():Tool(){}
+
+
+bool FilterEvents::Initialise(std::string configfile, DataModel &data){
+
+  /////////////////// Useful header ///////////////////////
+  if(configfile!="") m_variables.Initialise(configfile); // loading config file
+  //m_variables.Print();
+
+  m_data= &data; //assigning transient data pointer
+  /////////////////////////////////////////////////////////////////
+
+  m_variables.Get("verbosity",verbosity);
+  m_variables.Get("FilteredFilesBasename",FilteredFilesBasename);
+  m_variables.Get("SavePath",SavePath);
+
+  matched = 0;
+
+  FilteredEvents = new BoostStore(false,2);   
+
+  std::string EventSelectorCutConfiguration;
+  m_data->CStore.Get("CutConfiguration",EventSelectorCutConfiguration);
+  FilterName = EventSelectorCutConfiguration;
+
+  FilteredEvents->Header->Set("FilteredEvent",true);
+  FilteredEvents->Header->Set("FilterName",FilterName);
+
+  return true;
+}
+
+
+bool FilterEvents::Execute(){
+
+  std::map<std::string,bool> dataStreams;
+  int runNumber, eventNumber, localEventNumber, partNumber;
+  m_data->Stores["ANNIEEvent"]->Get("DataStreams",dataStreams);
+  m_data->Stores["ANNIEEvent"]->Get("RunNumber",runNumber);
+  m_data->Stores["ANNIEEvent"]->Get("EventNumber",eventNumber);
+  m_data->Stores["ANNIEEvent"]->Get("LocalEventNumber",localEventNumber);
+  m_data->Stores["ANNIEEvent"]->Get("PartNumber",partNumber);
+
+
+  bool pass_filter = false;
+  m_data->Stores.at("RecoEvent")->Get("EventCutStatus", pass_filter); 
+
+  if (pass_filter){
+    this->SetAndSaveEvent();
+    matched++;
+  }
+
+  return true;
+}
+
+
+bool FilterEvents::Finalise(){
+
+  FilteredEvents->Close();
+  FilteredEvents->Delete();
+  delete FilteredEvents;
+  std::cout<<"Number of matched run numbers, part numbers, and (local) event numbers: "<<matched<<std::endl;
+
+  return true;
+}
+
+void FilterEvents::SetAndSaveEvent(){
+
+  uint32_t RunNumber, SubrunNumber, EventNumber, TriggerWord;
+  size_t LocalEventNumber;
+  uint64_t RunStartTime, EventTimeTank, EventTimeLAPPD, CTCTimestamp, LAPPDOffset;
+  TimeClass EventTimeMRD;
+  int PartNumber, RunType, TriggerExtended;
+  std::map<std::string, bool> DataStreams;
+  std::string MRDTriggerType;
+  std::map<std::string, int> MRDLoopbackTDC;
+  TriggerClass TriggerData;
+  BeamStatus BeamStatus;
+  std::map<unsigned long, std::vector<Hit>>* TDCData;
+  std::map<unsigned long, std::vector<Hit>> *Hits, *AuxHits;
+  std::map<unsigned long, std::vector<int>> RawAcqSize;
+  std::map<unsigned long, std::vector<std::vector<ADCPulse>>> RecoADCData, RecoAuxADCData;
+  PsecData LAPPDData;
+
+  std::map<unsigned long, std::vector<Hit>> *NewTDCData = new std::map<unsigned long, std::vector<Hit>>;
+  std::map<unsigned long, std::vector<Hit>> *NewHits = new  std::map<unsigned long, std::vector<Hit>>;
+  std::map<unsigned long, std::vector<Hit>> *NewAuxHits = new std::map<unsigned long, std::vector<Hit>>;
+
+  // Get and Set all variables in the event to the new booststore
+  m_data->Stores["ANNIEEvent"]->Get("AuxHits",AuxHits);
+    for (auto&& entry : (*AuxHits)){
+        NewAuxHits->emplace(entry.first,entry.second);
+        }
+    FilteredEvents->Set("AuxHits",NewAuxHits,true);
+
+    m_data->Stores["ANNIEEvent"]->Get("BeamStatus",BeamStatus);  FilteredEvents->Set("BeamStatus",BeamStatus);
+    m_data->Stores["ANNIEEvent"]->Get("CTCTimestamp",CTCTimestamp);  FilteredEvents->Set("CTCTimestamp",CTCTimestamp);
+    m_data->Stores["ANNIEEvent"]->Get("DataStreams",DataStreams);  FilteredEvents->Set("DataStreams",DataStreams);
+    m_data->Stores["ANNIEEvent"]->Get("EventNumber",EventNumber);  FilteredEvents->Set("EventNumber",EventNumber);
+    m_data->Stores["ANNIEEvent"]->Get("EventTimeLAPPD",EventTimeLAPPD);  FilteredEvents->Set("EventTimeLAPPD",EventTimeLAPPD);
+    m_data->Stores["ANNIEEvent"]->Get("EventTimeMRD",EventTimeMRD);  FilteredEvents->Set("EventTimeMRD",EventTimeMRD);
+    m_data->Stores["ANNIEEvent"]->Get("EventTimeTank",EventTimeTank);  FilteredEvents->Set("EventTimeTank",EventTimeTank);
+    m_data->Stores["ANNIEEvent"]->Get("Hits",Hits);
+    for (auto&& entry : (*Hits)){
+        NewHits->emplace(entry.first,entry.second);
+        }
+    FilteredEvents->Set("Hits",NewHits,true);
+
+/*   
+    m_data->Stores["ANNIEEvent"]->Get("LAPPDData",LAPPDData);  FilteredEvents->Set("LAPPDData",LAPPDData);
+    m_data->Stores["ANNIEEvent"]->Get("LAPPDOffset",LAPPDOffset);  FilteredEvents->Set("LAPPDOffset",LAPPDOffset);
+*/
+    m_data->Stores["ANNIEEvent"]->Get("LocalEventNumber",LocalEventNumber);  FilteredEvents->Set("LocalEventNumber",LocalEventNumber);
+    m_data->Stores["ANNIEEvent"]->Get("MRDLoopbackTDC",MRDLoopbackTDC);  FilteredEvents->Set("MRDLoopbackTDC",MRDLoopbackTDC);
+    m_data->Stores["ANNIEEvent"]->Get("MRDTriggerType",MRDTriggerType);  FilteredEvents->Set("MRDTriggerType",MRDTriggerType);
+    m_data->Stores["ANNIEEvent"]->Get("PartNumber",PartNumber);  FilteredEvents->Set("PartNumber",PartNumber);
+    m_data->Stores["ANNIEEvent"]->Get("RawAcqSize",RawAcqSize);  FilteredEvents->Set("RawAcqSize",RawAcqSize);
+    m_data->Stores["ANNIEEvent"]->Get("RecoADCData",RecoADCData);  FilteredEvents->Set("RecoADCData",RecoADCData);
+    m_data->Stores["ANNIEEvent"]->Get("RecoAuxADCData",RecoAuxADCData);  FilteredEvents->Set("RecoAuxADCData",RecoAuxADCData);
+    m_data->Stores["ANNIEEvent"]->Get("RunNumber",RunNumber);  FilteredEvents->Set("RunNumber",RunNumber);
+    m_data->Stores["ANNIEEvent"]->Get("RunStartTime",RunStartTime);  FilteredEvents->Set("RunStartTime",RunStartTime);
+    m_data->Stores["ANNIEEvent"]->Get("RunType",RunType);    FilteredEvents->Set("RunType",RunType);
+    m_data->Stores["ANNIEEvent"]->Get("SubrunNumber",SubrunNumber);  FilteredEvents->Set("SubrunNumber",SubrunNumber);
+    m_data->Stores["ANNIEEvent"]->Get("TDCData",TDCData);
+    for (auto&& entry : (*TDCData)){
+        NewTDCData->emplace(entry.first,entry.second);
+        }
+
+    FilteredEvents->Set("TDCData",NewTDCData,true);
+    m_data->Stores["ANNIEEvent"]->Get("TriggerData",TriggerData);  FilteredEvents->Set("TriggerData",TriggerData);
+    m_data->Stores["ANNIEEvent"]->Get("TriggerExtended",TriggerExtended);  FilteredEvents->Set("TriggerExtended",TriggerExtended);
+    m_data->Stores["ANNIEEvent"]->Get("TriggerWord",TriggerWord);  FilteredEvents->Set("TriggerWord",TriggerWord);
+
+    std::string Filename = SavePath + "/" + FilteredFilesBasename + "_" + FilterName;
+
+    if (verbosity>0) std::cout<<"Filename is "<<Filename<<std::endl;
+    FilteredEvents->Save(Filename);
+    FilteredEvents->Delete();
+
+}

UserTools/FilterEvents/FilterEvents.h

@@ -0,0 +1,52 @@
+#ifndef FilterEvents_H
+#define FilterEvents_H
+
+#include <string>
+#include <iostream>
+
+#include "Tool.h"
+
+#include "ADCPulse.h"
+#include "PsecData.h"
+#include "Hit.h"
+
+/**
+ * \class FilterEvents
+ *
+ * This is a blank template for a Tool used by the script to generate a new custom tool. Please fill out the description and author information.
+*
+* $Author: B.Richards $
+* $Date: 2019/05/28 10:44:00 $
+* Contact: [email protected]
+*/
+class FilterEvents: public Tool {
+
+
+ public:
+
+  FilterEvents(); ///< Simple constructor
+  bool Initialise(std::string configfile,DataModel &data); ///< Initialise Function for setting up Tool resources. @param configfile The path and name of the dynamic configuration file to read in. @param data A reference to the transient data class used to pass information between Tools.
+  bool Execute(); ///< Execute function used to perform Tool purpose.
+  bool Finalise(); ///< Finalise function used to clean up resources.
+
+  void SetAndSaveEvent();
+
+ private:
+
+  std::string FilterName;
+  std::string FilteredFilesBasename;
+  std::string SavePath;
+  int verbosity;
+  BoostStore* FilteredEvents = nullptr;
+  int matched;
+
+  int v_error = 0;
+  int v_warning = 1;
+  int v_message = 2;
+  int v_debug = 3;
+  int vv_debug = 4;
+
+};
+
+
+#endif

UserTools/FilterEvents/README.md

@@ -0,0 +1,26 @@
+# FilterEvents
+
+The `FilterEvents` tool reads in ANNIEEvent BoostStore files and saves them to a Filtered output BoostStore files in case selected cuts were passed for a given event. The cuts themselves are defined in the `EventSelector` tool and official filters should be defined in the `configfiles/Filters/` directory.
+
+## Data
+
+The `FilterEvents` tool reads in the complete `ANNIEEvent` BoostStore file and saves all the objects to a new, filtered output BoostStore file. 
+
+* Input data: `ANNIEEvent` BoostStore (typically loaded with the `LoadANNIEEvent` tool)
+* Output data: `ANNIEEvent` BoostStore, saved to a new file. Only events which passed the cuts specified in the `EventSelector` tool are saved to this output BoostStore file
+
+Within the output BoostStore file, a separate flag indicates that the events have been filtered. This flag is a Boolean and is called "FilteredEvent". In addition, the FilterName is also saved to the output BoostStore:
+
+* `FilteredEvent (bool)`: indicates whether the object is part of a filtered output BoostStore file
+* `FilterName (string)`: The name of the applied filter in the `EventSelector` tool.
+
+## Configuration
+
+`FilterEvents` uses the following configuration variables:
+
+```
+verbosity 2                             #how much logging output do you want?
+FilteredFilesBasename FilteredEvents    #base name of the BoostStore output files
+                                        #the name of the Filter will be appended automatically to this basename
+SavePath /path/to/output/               #where should the output files be stored
+```

UserTools/LoadANNIEEvent/LoadANNIEEvent.cpp

@@ -271,13 +271,15 @@ bool LoadANNIEEvent::Execute() {
   if ((int)current_entry_ != offset_evnum) m_data->Stores["ANNIEEvent"]->Delete();	//ensures that we can access pointers without problems
 
   m_data->Stores["ANNIEEvent"]->GetEntry(current_entry_);  
-  if (!m_data->Stores["ANNIEEvent"]->Has("LocalEventNumber")){ m_data->Stores["ANNIEEvent"]->Set("LocalEventNumber",current_entry_); }
+  bool has_local = (m_data->Stores["ANNIEEvent"]->Has("LocalEventNumber"));
+  bool has_global = (m_data->Stores["ANNIEEvent"]->Has("EventNumber"));
+  if (!has_local){ m_data->Stores["ANNIEEvent"]->Set("LocalEventNumber",current_entry_);}
   ++current_entry_;
 
-  if (global_evnr && !m_data->Stores["ANNIEEvent"]->Has("LocalEventNumber")){ m_data->Stores["ANNIEEvent"]->Set("EventNumber",global_ev); }
-  global_ev++; 
 
 
+  if (global_evnr && !has_local){ m_data->Stores["ANNIEEvent"]->Set("EventNumber",global_ev); }
+  global_ev++; 
 
   if ( current_entry_ >= total_entries_in_file_ ) {
     ++current_file_;

UserTools/Unity.h

@@ -164,4 +164,5 @@
 #include "PythonScript.h"
 #include "ReweightEventsGenie.h"
 #include "FilterLAPPDEvents.h"
+#include "FilterEvents.h"
 #include "Stage1DataBuilder.h"

configfiles/BeamClusterAnalysis/PhaseIITreeMakerConfig

@@ -13,4 +13,5 @@ MCTruth_fill 0
 Reco_fill 0
 RecoDebug_fill 0
 muonTruthRecoDiff_fill 0
-
+IsData 1
+HasGenie 0

configfiles/FilterEvents/ClusterClassifiersConfig

@@ -0,0 +1 @@
+verbosity 0

configfiles/FilterEvents/ClusterFinderConfig

@@ -0,0 +1,12 @@
+# ClusterFinder Config File
+
+verbosity 0
+HitStore Hits #Either MCHits or Hits (accessed in ANNIEEvent store)
+OutputFile BeamRun_ClusterFinder_DefaultOutput #Output root prefix name for the current run
+ClusterFindingWindow 100 # in ns, size of the window used to "clusterize"
+AcqTimeWindow 70000 # in ns, size of the acquisition window
+ClusterIntegrationWindow 100 # in ns, all hits with +/- 1/2 of this window are considered in the cluster
+MinHitsPerCluster 5 # group of hits are considered clusters above this amount of hits
+end_of_window_time_cut 0.95 # from o to 1, length of the window you want to loop over with respect to acq. window (1 for full window, 0.95 for 95% from the start)
+Plots2D 0 #Draw 2D charge-vs-time plots?
+ChankeyToPMTIDMap ./configfiles/EventDisplay/Data-RecoEvent/Chankey_WCSimID.dat

configfiles/FilterEvents/FilterEventsConfig

@@ -0,0 +1,5 @@
+# FilterEvents config file
+
+verbosity 2
+FilteredFilesBasename FilteredEvents_R2630
+SavePath .

configfiles/FilterEvents/FindMrdTracksConfig

@@ -0,0 +1,12 @@
+# FindMrdTracks Config File
+# all variables retrieved with m_variables.Get() must be defined here!
+
+verbosity 0
+IsData 1
+OutputDirectory .
+OutputFile STEC_MRDTracks_cluster40ns
+DrawTruthTracks 0 # whether to add MC Truth track info for drawing in MrdPaddlePlot Tool
+                  ## note you need to run that tool to actually view the tracks!
+WriteTracksToFile 0	# should the track information be written to a ROOT-file?
+SelectTriggerType 0	#should the loaded data be filtered by trigger type?
+TriggerType Beam	#options: Cosmic, Beam, No Loopback

configfiles/FilterEvents/LoadANNIEEventConfig

@@ -0,0 +1,4 @@
+verbose 0
+FileForListOfInputs ./configfiles/FilterEvents/my_inputs.txt
+EventOffset 0
+GlobalEvNr 1

configfiles/FilterEvents/README.md

@@ -0,0 +1,38 @@
+# FilterEvents toolchain
+
+***********************
+# Description
+**********************
+
+The `FilterEvents` toolchain uses cuts from the `EventSelector` tool to filter out events passing the selected cuts. The filtered events are then saved to a new dedicated output BoostStore that only contains the filtered events. The output BoostStore file will contain all the variables from the `ANNIEEvent` BoostStore and can be used in a toolchain like a regular processed data file.
+
+Filters that define the specific event selection cuts should be defined in the `configfiles/Filters/` directory. Currently the following filters are defined:
+
+* NeutrinoCandidate
+* DirtCandidate
+* Throughgoing
+* BeamTrigger
+* CosmicTrigger
+
+Each filter is defined by a unique name, which is passed on to the `FilterEvents` tool via the common `CStore`. The variable is then used as a part of the output BoostStore file. This guarantees that the output BoostStore will contain the name of the official Filter and it is easily recognizable which filter was used for a specific file.
+
+************************
+# ToolsConfig
+************************
+
+The following tools are in the `FilterEvents` toolchain:
+
+```
+myLoadANNIEEvent LoadANNIEEvent ./configfiles/FilterEvents/LoadANNIEEventConfig
+myLoadGeometry LoadGeometry ./configfiles/LoadGeometry/LoadGeometryConfig
+myTimeClustering TimeClustering configfiles/FilterEvents/TimeClusteringConfig
+myFindMrdTracks FindMrdTracks configfiles/FilterEvents/FindMrdTracksConfig
+myClusterFinder ClusterFinder ./configfiles/FilterEvents/ClusterFinderConfig
+myClusterClassifiers ClusterClassifiers ./configfiles/FilterEvents/ClusterClassifiersConfig
+myEventSelector EventSelector ./configfiles/Filters/NeutrinoCandidateFilter
+myFilterEvents FilterEvents ./configfiles/FilterEvents/FilterEventsConfig
+```
+
+The `EventSelector` row can be adjusted if a different filter should be used.
+
+New filters can be defined in the `configfiles/Filters` directory. However, old filters should not be changed in any way to preserve the uniqueness of filters. If a filter needs to be updated, one should rather define a new filter with a unique name, e.g. `NeutrinoCandidate_v2` for updated neutrino candidate cuts.

configfiles/FilterEvents/TimeClusteringConfig

@@ -0,0 +1,13 @@
+#TimeClustering config file
+
+verbosity 0
+MinDigitsForTrack 3
+MaxMrdSubEventDuration 30
+MinSubeventTimeSep 30
+MakeMrdDigitTimePlot 0
+LaunchTApplication 0
+IsData 1
+#OutputROOTFile TimeClustering_MRDTest28_cluster40ns
+OutputROOTFile STEC_TimeClusteringOut
+MapChankey_WCSimID ./configfiles/SimpleTankEnergyCalibrator/MRD_Chankey_WCSimID.dat
+