All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions | Private Attributes | List of all members
caldata::CalWire Class Reference
Inheritance diagram for caldata::CalWire:

Public Member Functions

 CalWire (fhicl::ParameterSet const &pset)
 
void produce (art::Event &evt)
 
void beginJob ()
 

Private Attributes

std::string fResponseFile
 
int fExpEndBins
 number of end bins to consider for tail fit More...
 
int fPostsample
 number of postsample bins More...
 
std::string fDigitModuleLabel
 module that made digits More...
 
std::vector< std::vector
< TComplex > > 
fKernelR
 
std::vector< std::vector
< TComplex > > 
fKernelS
 
std::vector< double > fDecayConstsR
 
std::vector< double > fDecayConstsS
 
std::vector< int > fKernMapR
 
std::vector< int > fKernMapS
 

Detailed Description

Definition at line 42 of file CalWire_module.cc.

Constructor & Destructor Documentation

caldata::CalWire::CalWire ( fhicl::ParameterSet const &  pset)
explicit

Definition at line 83 of file CalWire_module.cc.

84  : EDProducer{pset}
85  {
86  fDigitModuleLabel = pset.get< std::string >("DigitModuleLabel", "daq");
87  cet::search_path sp("FW_SEARCH_PATH");
88  sp.find_file(pset.get<std::string>("ResponseFile"), fResponseFile);
89  fExpEndBins = pset.get< int > ("ExponentialEndBins");
90  fPostsample = pset.get< int > ("PostsampleBins");
91 
92  produces< std::vector<recob::Wire> >();
93  produces<art::Assns<raw::RawDigit, recob::Wire>>();
94 
95  }
std::string fResponseFile
std::string fDigitModuleLabel
module that made digits
int fExpEndBins
number of end bins to consider for tail fit
int fPostsample
number of postsample bins

Member Function Documentation

void caldata::CalWire::beginJob ( )

Definition at line 98 of file CalWire_module.cc.

99  {
100 
101  MF_LOG_DEBUG("CalWire") << "CalWire_plugin: Opening Electronics Response File: "
102  << fResponseFile.c_str();
103 
104  TFile f(fResponseFile.c_str());
105  if( f.IsZombie() )
106  mf::LogWarning("CalWire") << "Cannot open response file "
107  << fResponseFile.c_str();
108 
109  TH2D *respRe = dynamic_cast<TH2D*>(f.Get("real/RespRe") );
110  TH2D *respIm = dynamic_cast<TH2D*>(f.Get("real/RespIm") );
111  TH1D *decayHist = dynamic_cast<TH1D*>(f.Get("real/decayHist"));
112  unsigned int wires = decayHist->GetNbinsX();
113  unsigned int bins = respRe->GetYaxis()->GetNbins();
114  unsigned int bin = 0;
115  unsigned int wire = 0;
116  fDecayConstsR.resize(wires);
117  fKernMapR.resize(wires);
118  fKernelR.resize(respRe->GetXaxis()->GetNbins());
119  const TArrayD *edges = respRe->GetXaxis()->GetXbins();
120  for(int i = 0; i < respRe->GetXaxis()->GetNbins(); ++i) {
121  fKernelR[i].resize(bins);
122  for(bin = 0; bin < bins; ++bin) {
123 
124  const TComplex a(respRe->GetBinContent(i+1,bin+1),
125  respIm->GetBinContent(i+1,bin+1));
126  fKernelR[i][bin]=a;
127  }
128  for(; wire < (*edges)[i+1]; ++wire) {
129  fKernMapR[wire]=i;
130  fDecayConstsR[wire]=decayHist->GetBinContent(wire+1);
131  }
132  }
133  respRe = dynamic_cast<TH2D*>(f.Get("sim/RespRe") );
134  respIm = dynamic_cast<TH2D*>(f.Get("sim/RespIm") );
135  decayHist = dynamic_cast<TH1D*>(f.Get("sim/decayHist"));
136  wires = decayHist->GetNbinsX();
137  bins = respRe->GetYaxis()->GetNbins();
138  fDecayConstsS.resize(wires);
139  fKernMapS.resize(wires);
140  fKernelS.resize(respRe->GetXaxis()->GetNbins());
141  const TArrayD *edges1 = respRe->GetXaxis()->GetXbins();
142  wire =0;
143  for(int i = 0; i < respRe->GetXaxis()->GetNbins(); ++i) {
144  fKernelS[i].resize(bins);
145  for(bin = 0; bin < bins; ++bin) {
146  const TComplex b(respRe->GetBinContent(i+1,bin+1),
147  respIm->GetBinContent(i+1,bin+1));
148  fKernelS[i][bin]=b;
149  }
150  for(; wire < (*edges1)[i+1]; ++wire) {
151  fKernMapS[wire]=i;
152  fDecayConstsS[wire]=decayHist->GetBinContent(wire+1);
153  }
154  }
155 
156  f.Close();
157  }
std::string fResponseFile
std::vector< int > fKernMapR
std::vector< double > fDecayConstsR
constexpr details::BinObj< T > bin(T value)
Returns a wrapper to print the specified data in binary format.
process_name gaushit a
std::vector< double > fDecayConstsS
std::vector< std::vector< TComplex > > fKernelR
std::vector< std::vector< TComplex > > fKernelS
std::vector< int > fKernMapS
void caldata::CalWire::produce ( art::Event &  evt)

Definition at line 160 of file CalWire_module.cc.

161  {
162 
163 
164  // get the geometry
165  art::ServiceHandle<geo::Geometry const> geom;
166 
167  std::vector<double> decayConsts;
168  std::vector<int> kernMap;
169  std::vector<std::vector<TComplex> > kernel;
170  //Put correct response functions and decay constants in place
171  if(evt.isRealData()) {
172  decayConsts=fDecayConstsR;
173  kernMap=fKernMapR;
174  kernel=fKernelR;
175  }
176  else {
177  decayConsts=fDecayConstsS;
178  kernMap=fKernMapS;
179  kernel=fKernelS;
180  }
181 
182  // get the FFT service to have access to the FFT size
183  art::ServiceHandle<util::LArFFT> fFFT;
184 
185  // make a collection of Wires
186  std::unique_ptr<std::vector<recob::Wire> > wirecol(new std::vector<recob::Wire>);
187  // ... and an association set
188  std::unique_ptr<art::Assns<raw::RawDigit,recob::Wire> > WireDigitAssn
189  (new art::Assns<raw::RawDigit,recob::Wire>);
190 
191  // Read in the digit List object(s).
192  art::Handle< std::vector<raw::RawDigit> > digitVecHandle;
193  evt.getByLabel(fDigitModuleLabel, digitVecHandle);
194 
195  if (!digitVecHandle->size()) return;
196  mf::LogInfo("CalWire") << "CalWire:: digitVecHandle size is " << digitVecHandle->size();
197 
198  // Use the handle to get a particular (0th) element of collection.
199  art::Ptr<raw::RawDigit> digitVec0(digitVecHandle, 0);
200 
201  unsigned int dataSize = digitVec0->Samples(); //size of raw data vectors
202 
203  int transformSize = fFFT->FFTSize();
204  raw::ChannelID_t channel(raw::InvalidChannelID); // channel number
205  unsigned int bin(0); // time bin loop variable
206 
207  double decayConst = 0.; // exponential decay constant of electronics shaping
208  double fitAmplitude = 0.; //This is the seed value for the amplitude in the exponential tail fit
209  std::vector<float> holder; // holds signal data
210  std::vector<short> rawadc(transformSize); // vector holding uncompressed adc values
211  std::vector<TComplex> freqHolder(transformSize+1); // temporary frequency data
212 
213  // loop over all wires
214  for(unsigned int rdIter = 0; rdIter < digitVecHandle->size(); ++rdIter){ // ++ move
215  holder.clear();
216 
217  art::Ptr<raw::RawDigit> digitVec(digitVecHandle, rdIter);
218  channel = digitVec->Channel();
219 
220  holder.resize(transformSize);
221 
222  // uncompress the data
223  raw::Uncompress(digitVec->ADCs(), rawadc, digitVec->Compression());
224 
225  for(bin = 0; bin < dataSize; ++bin)
226  holder[bin]=(rawadc[bin]-digitVec->GetPedestal());
227  // fExpEndBins only nonzero for detectors needing exponential tail fitting
228  if(fExpEndBins && std::abs(decayConsts[channel]) > 0.0){
229 
230  TH1D expTailData("expTailData","Tail data for fit",
231  fExpEndBins,dataSize-fExpEndBins,dataSize);
232  TF1 expFit("expFit","[0]*exp([1]*x)");
233 
234  for(bin = 0; bin < (unsigned int)fExpEndBins; ++bin)
235  expTailData.Fill(dataSize-fExpEndBins+bin,holder[dataSize-fExpEndBins+bin]);
236  decayConst = decayConsts[channel];
237  fitAmplitude = holder[dataSize-fExpEndBins]/exp(decayConst*(dataSize-fExpEndBins));
238  expFit.FixParameter(1,decayConst);
239  expFit.SetParameter(0,fitAmplitude);
240  expTailData.Fit(&expFit,"QWN","",dataSize-fExpEndBins,dataSize);
241  expFit.SetRange(dataSize,transformSize);
242  for(bin = 0; bin < dataSize; ++bin)
243  holder[dataSize+bin]= expFit.Eval(bin+dataSize);
244  }
245  // This is actually deconvolution, by way of convolution with the inverted
246  // kernel. This code assumes the response function has already been
247  // been transformed and inverted. This way a complex multiplication, rather
248  // than a complex division is performed saving 2 multiplications and
249  // 2 divsions
250 
251  // the example below is for MicroBooNE, experiments should
252  // adapt as appropriate
253 
254  // Figure out which kernel to use (0=induction, 1=collection).
255  geo::SigType_t sigtype = geom->SignalType(channel);
256  size_t k;
257  if(sigtype == geo::kInduction)
258  k = 0;
259  else if(sigtype == geo::kCollection)
260  k = 1;
261  else
262  throw cet::exception("CalWire") << "Bad signal type = " << sigtype << "\n";
263  if (k >= kernel.size())
264  throw cet::exception("CalWire") << "kernel size < " << k << "!\n";
265 
266  fFFT->Convolute(holder,kernel[k]);
267 
268  holder.resize(dataSize,1e-5);
269  //This restores the DC component to signal removed by the deconvolution.
270  if(fPostsample) {
271  double average=0.0;
272  for(bin=0; bin < (unsigned int)fPostsample; ++bin)
273  average+=holder[holder.size()-1-bin]/(double)fPostsample;
274  for(bin = 0; bin < holder.size(); ++bin) holder[bin]-=average;
275  }
276  wirecol->push_back(recob::WireCreator(holder,*digitVec).move());
277  // add an association between the last object in wirecol
278  // (that we just inserted) and digitVec
279  if (!util::CreateAssn(*this, evt, *wirecol, digitVec, *WireDigitAssn)) {
280  throw art::Exception(art::errors::ProductRegistrationFailure)
281  << "Can't associate wire #" << (wirecol->size() - 1)
282  << " with raw digit #" << digitVec.key();
283  } // if failed to add association
284  } // for raw digits
285 
286  if(wirecol->size() == 0)
287  mf::LogWarning("CalWire") << "No wires made for this event.";
288 
289  evt.put(std::move(wirecol));
290  evt.put(std::move(WireDigitAssn));
291 
292  return;
293  }
std::vector< int > fKernMapR
std::vector< double > fDecayConstsR
Class managing the creation of a new recob::Wire object.
Definition: WireCreator.h:53
constexpr details::BinObj< T > bin(T value)
Returns a wrapper to print the specified data in binary format.
std::vector< double > fDecayConstsS
T abs(T value)
constexpr ChannelID_t InvalidChannelID
ID of an invalid channel.
Definition: RawTypes.h:32
Signal from induction planes.
Definition: geo_types.h:145
enum geo::_plane_sigtype SigType_t
std::vector< std::vector< TComplex > > fKernelR
bool CreateAssn(art::Event &evt, std::vector< T > const &a, art::Ptr< U > const &b, art::Assns< U, T > &assn, std::string a_instance, size_t index=UINT_MAX)
Creates a single one-to-one association.
std::string fDigitModuleLabel
module that made digits
std::vector< std::vector< TComplex > > fKernelS
int fExpEndBins
number of end bins to consider for tail fit
do i e
int fPostsample
number of postsample bins
TCEvent evt
Definition: DataStructs.cxx:8
pdgs k
Definition: selectors.fcl:22
std::vector< int > fKernMapS
unsigned int ChannelID_t
Type representing the ID of a readout channel.
Definition: RawTypes.h:28
void Uncompress(const std::vector< short > &adc, std::vector< short > &uncompressed, raw::Compress_t compress)
Uncompresses a raw data buffer.
Definition: raw.cxx:776
Signal from collection planes.
Definition: geo_types.h:146

Member Data Documentation

std::vector<double> caldata::CalWire::fDecayConstsR
private

vector holding RC decay constants

Definition at line 67 of file CalWire_module.cc.

std::vector<double> caldata::CalWire::fDecayConstsS
private

vector holding RC decay constants

Definition at line 69 of file CalWire_module.cc.

std::string caldata::CalWire::fDigitModuleLabel
private

module that made digits

Definition at line 61 of file CalWire_module.cc.

int caldata::CalWire::fExpEndBins
private

number of end bins to consider for tail fit

Definition at line 59 of file CalWire_module.cc.

std::vector<std::vector<TComplex> > caldata::CalWire::fKernelR
private

holds transformed induction response function

Definition at line 63 of file CalWire_module.cc.

std::vector<std::vector<TComplex> > caldata::CalWire::fKernelS
private

holds transformed induction response function

Definition at line 65 of file CalWire_module.cc.

std::vector<int> caldata::CalWire::fKernMapR
private

map telling which channels have which response functions

Definition at line 71 of file CalWire_module.cc.

std::vector<int> caldata::CalWire::fKernMapS
private

map telling which channels have which response functions

Definition at line 73 of file CalWire_module.cc.

int caldata::CalWire::fPostsample
private

number of postsample bins

Definition at line 60 of file CalWire_module.cc.

std::string caldata::CalWire::fResponseFile
private

response file containing transformed shape histograms and decay constants

Definition at line 55 of file CalWire_module.cc.


The documentation for this class was generated from the following file: