MergeDaqCollections_module.cc
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////
2 // \file MergeDaqCollections_module.cxx
3 // \brief Merge collections.
4 // At this step collection from primary and secondary files are merged under the same label. You can set RawDigit::SetMC(true or false). Which will let the calibrator use mc or real calibration constants when you analyze the merged file.
5 // I learned from MergeCollections_module.cc
6 // \author Enhao Song es4re@virginia.edu
7 ///////////////////////////////////////////////////////////////////////////
8 
9 #include <math.h>
10 #include <string>
11 #include <algorithm>
12 #include <limits>
13 
17 #include "fhiclcpp/ParameterSet.h"
21 
22 #include "DAQChannelMap/DAQChannelMap.h"
23 #include "RawData/RawDigit.h"
24 #include "RawData/RawTrigger.h"
25 #include "RawData/DAQHeader.h"
26 #include "RecoBase/CellHit.h"
27 #include "Calibrator/Calibrator.h"
28 #include "MCCheater/BackTracker.h"
29 
30 #include "Utilities/AssociationUtil.h"
31 
32 ///group hits according to time and space
34 
36  public:
37  explicit MergeDaqCollections(fhicl::ParameterSet const &pset);
38  virtual ~MergeDaqCollections();
39 
40  void produce(art::Event& evt);
41 
42  protected:
43  // internal methods
44 
45  // configuration settings
46  std::string fFirstDaqCollection; /// Label on the main collection
47  std::string fSecondCollection; /// Label on the main collection
48 
49  bool fFirstCollectionIsMC; /// First collection is MC
50  bool fSecondCollectionIsMC; /// Second collection is MC
51  bool fUseSecondaryTrigger; /// Want trigger from 2nd collection rather than 1st
52 
53  bool fRemoveFirstCollectionNoise; /// Don't mix in noise hits from first collection, if it is MC.
54  bool fRemoveSecondCollectionNoise; /// Don't mix in noise hits from second collection, if it is MC.
55  bool fOverlayWithTDCCut; /// Only overlay secondary hits within a TDC range of the main event hits. Since we are using the secondary hits as background noise hits. We only need them in the time window close to the signal hits. We can still use the cosmic files before overlay for as a larger data set for background study. This way we can save a lot of disk space for the overlaid files and computing time in the following analysis chain.
57 
58  std::string fInstanceLabel; /// Instance label
59  };
60 
61 
62  //----------------------------------------------------------------------
64  fFirstDaqCollection (pset.get< std::string >("FirstDaqCollection")),
65  fSecondCollection (pset.get< std::string >("SecondCollection")),
66  fFirstCollectionIsMC (pset.get< bool >("FirstCollectionIsMC")),
67  fSecondCollectionIsMC (pset.get< bool >("SecondCollectionIsMC")),
68  fUseSecondaryTrigger (pset.get< bool >("UseSecondaryTrigger")),
69  fRemoveFirstCollectionNoise (pset.get< bool >("RemoveFirstCollectionNoise")),
70  fRemoveSecondCollectionNoise(pset.get< bool >("RemoveSecondCollectionNoise")),
71  fOverlayWithTDCCut(pset.get< bool >("OverlayWithTDCCut")),
72  fTDCExtendValue(pset.get< int >("TDCExtendValue")),
73  fInstanceLabel (pset.get< std::string >("InstanceLabel"))
74  {
75  produces< std::vector<rawdata::RawDigit> >(fInstanceLabel);
76  produces< std::vector<rawdata::RawTrigger> >(fInstanceLabel);
77  }
78 
79  //----------------------------------------------------------------------
80  MergeDaqCollections::~MergeDaqCollections()
81  {
82  }
83 
84  //----------------------------------------------------------------------
85  void MergeDaqCollections::produce(art::Event& evt)
86  {
89 
90  // BackTracker has not yet seen MC mixed from secondary file
91  if(fSecondCollectionIsMC) bt->Rebuild(evt);
92 
93  std::unique_ptr< std::vector<rawdata::RawDigit> > digitcol (new std::vector<rawdata::RawDigit>);
94  std::unique_ptr< std::vector<rawdata::RawTrigger> > triggercol (new std::vector<rawdata::RawTrigger>);
95 
97  evt.getByLabel(fFirstDaqCollection, digit1list);
98 
100  evt.getByLabel(fSecondCollection, digit2list);
101 
102 
103  // First Collection Digits
104  auto min_TDC = std::numeric_limits<int32_t>::max();
105  auto max_TDC = std::numeric_limits<int32_t>::min();
106  for(unsigned int i = 0; i < digit1list->size(); ++i){
107  rawdata::RawDigit digit ((*digit1list)[i]);
108  digit.SetMC(fFirstCollectionIsMC);
109  if(fFirstCollectionIsMC && fRemoveFirstCollectionNoise){
110  // Make cell hit so backtracker can
111  // attempt to match to particle
112  if(!bt->IsNoise(cal->MakeCellHit(&digit)))
113  digitcol->push_back(digit);
114  }
115  else
116  digitcol->push_back(digit);
117  }
118  //find max and min TDC of the digitcol
119  if (fOverlayWithTDCCut){
120  for(auto digit: *digitcol){
121  min_TDC = std::min(min_TDC, digit.TDC());
122  max_TDC = std::max(max_TDC, digit.TDC());
123  }
124  }
125 
126  //extend time window on each end with fTDCExtendValue
127  min_TDC = std::max(min_TDC-fTDCExtendValue, -3200); // data always has TDC > -3200 TDC which is 50 us, so we only overlay data beyond this time point.
128  max_TDC = max_TDC+fTDCExtendValue;
129 
130  // Second Collection Digits
131  for(unsigned int i = 0; i < digit2list->size(); ++i){
132  rawdata::RawDigit digit ((*digit2list)[i]);
133  if (fOverlayWithTDCCut){
134  if (digit.TDC()>max_TDC || digit.TDC()<min_TDC)
135  continue;
136  }
137  digit.SetMC(fSecondCollectionIsMC);
138  if(fSecondCollectionIsMC && fRemoveSecondCollectionNoise){
139  // Make cell hit so backtracker can
140  // attempt to match to particle
141  if(!bt->IsNoise(cal->MakeCellHit(&digit)))
142  digitcol->push_back(digit);
143  }
144  else
145  digitcol->push_back(digit);
146  }
147 
148 
149  // It can't be assumed which trigger is appropriate
151 
152  if(fUseSecondaryTrigger){
153  evt.getByLabel(fSecondCollection, triggerlist);
154 
155  for (unsigned int i = 0; i < triggerlist->size(); ++i) {
156  rawdata::RawTrigger trigger ((*triggerlist)[i]);
157  triggercol->push_back(trigger);
158  }
159  } else {
160  evt.getByLabel(fFirstDaqCollection, triggerlist);
161 
162  for (unsigned int i = 0; i < triggerlist->size(); ++i) {
163  rawdata::RawTrigger trigger ((*triggerlist)[i]);
164  triggercol->push_back(trigger);
165  }
166  }
167 
168  // put the collections in the event
169  evt.put(std::move(digitcol), fInstanceLabel);
170  evt.put(std::move(triggercol), fInstanceLabel);
171 
172  }
173 
175 } // end namespace MergeDaqCollections
176 /////////////////////////////////////////////////////////////////////////
T max(const caf::Proxy< T > &a, T b)
back track the reconstruction to the simulation
int32_t TDC() const
The time of the last baseline sample.
Definition: RawDigit.h:94
bool fOverlayWithTDCCut
Don&#39;t mix in noise hits from second collection, if it is MC.
DEFINE_ART_MODULE(TestTMapFile)
rb::CellHit MakeCellHit(const rawdata::RawDigit *rawdigit)
bool IsNoise(const art::Ptr< rb::CellHit > &hit) const
Is this hit not associated with any particles?
ProductID put(std::unique_ptr< PROD > &&product)
Definition: Event.h:102
bool fRemoveSecondCollectionNoise
Don&#39;t mix in noise hits from first collection, if it is MC.
bool fFirstCollectionIsMC
Label on the main collection.
int evt
std::string fSecondCollection
Label on the main collection.
bool fRemoveFirstCollectionNoise
Want trigger from 2nd collection rather than 1st.
void Rebuild(const art::Event &evt)
void SetMC(bool isMC=true)
Definition: RawDigit.h:106
static float min(const float a, const float b, const float c)
Definition: absgeo.cxx:45
bool getByLabel(std::string const &label, std::string const &productInstanceName, Handle< PROD > &result) const
Definition: DataViewImpl.h:344
int fTDCExtendValue
Only overlay secondary hits within a TDC range of the main event hits. Since we are using the seconda...
group hits according to time and space
T min(const caf::Proxy< T > &a, T b)
T max(sqlite3 *const db, std::string const &table_name, std::string const &column_name)
Definition: statistics.h:68
enum BeamMode string