util::SignalShapingServiceDUNE35t Class Reference

#include <SignalShapingServiceDUNE35t.h>

Public Member Functions
	SignalShapingServiceDUNE35t (const fhicl::ParameterSet &pset, art::ActivityRegistry &reg)
	~SignalShapingServiceDUNE35t ()
void	reconfigure (const fhicl::ParameterSet &pset)
std::vector< DoubleVec >	GetNoiseFactVec ()
double	GetASICGain (unsigned int const channel) const
double	GetShapingTime (unsigned int const channel) const
double	GetRawNoise (unsigned int const channel) const
double	GetDeconNoise (unsigned int const channel) const
const util::SignalShaping &	SignalShaping (unsigned int channel) const
int	FieldResponseTOffset (detinfo::DetectorClocksData const &clockData, unsigned int const channel) const
template<class T >
void	Convolute (detinfo::DetectorClocksData const &clockData, unsigned int channel, std::vector< T > &func) const
template<class T >
void	Deconvolute (detinfo::DetectorClocksData const &clockData, unsigned int channel, std::vector< T > &func) const
double	GetDeconNorm ()

Private Member Functions
void	init () const
void	init ()
void	SetFieldResponse (detinfo::DetectorClocksData const &clockData, detinfo::DetectorPropertiesData const &detProp)
void	SetElectResponse (double shapingtime, double gain)
void	SetFilters (detinfo::DetectorClocksData const &clockData)
void	SetResponseSampling (detinfo::DetectorClocksData const &clockData)

Private Attributes
bool	fInit
	Initialization flag. More...
int	fNFieldBins
	number of bins for field response More...
double	fCol3DCorrection
double	fInd3DCorrection
double	fColFieldRespAmp
	amplitude of response to field More...
double	fIndUFieldRespAmp
	amplitude of response to field More...
double	fIndVFieldRespAmp
	amplitude of response to field More...
std::vector< double >	fFieldResponseTOffset
	Time offset for field response in ns. More...
std::vector< double >	fCalibResponseTOffset
double	fInputFieldRespSamplingPeriod
	Sampling period in the input field response. More...
double	fDeconNorm
double	fADCPerPCAtLowestASICGain
	Pulse amplitude gain for a 1 pc charge impulse after convoluting it the with field and electronics response with the lowest ASIC gain setting of 4.7 mV/fC. More...
std::vector< DoubleVec >	fNoiseFactVec
std::vector< double >	fASICGainInMVPerFC
std::vector< double >	fShapeTimeConst
	time constants for exponential shaping More...
TF1 *	fColFilterFunc
	Parameterized collection filter function. More...
TF1 *	fIndUFilterFunc
	Parameterized induction filter function. More...
TF1 *	fIndVFilterFunc
	Parameterized induction filter function. More...
bool	fUseFunctionFieldShape
	Flag that allows to use a parameterized field response instead of the hardcoded version. More...
bool	fUseHistogramFieldShape
	Flag that turns on field response shapes from histograms. More...
bool	fGetFilterFromHisto
	Flag that allows to use a filter function from a histogram instead of the functional dependency. More...
TF1 *	fColFieldFunc
	Parameterized collection field shape function. More...
TF1 *	fIndUFieldFunc
	Parameterized induction field shape function. More...
TF1 *	fIndVFieldFunc
	Parameterized induction field shape function. More...
TH1F *	fFieldResponseHist [3]
	Histogram used to hold the field response, hardcoded for the time being. More...
TH1D *	fFilterHist [3]
	Histogram used to hold the collection filter, hardcoded for the time being. More...
util::SignalShaping	fColSignalShaping
util::SignalShaping	fIndUSignalShaping
util::SignalShaping	fIndVSignalShaping
std::vector< double >	fColFieldResponse
std::vector< double >	fIndUFieldResponse
std::vector< double >	fIndVFieldResponse
std::vector< double >	fElectResponse
std::vector< TComplex >	fColFilter
std::vector< TComplex >	fIndUFilter
std::vector< TComplex >	fIndVFilter
bool	fDebugFieldShape
	Save field response shape to a text file. More...

Detailed Description

Definition at line 48 of file SignalShapingServiceDUNE35t.h.

Constructor & Destructor Documentation

util::SignalShapingServiceDUNE35t::SignalShapingServiceDUNE35t	(	const fhicl::ParameterSet &	pset,
		art::ActivityRegistry &	reg
	)

Definition at line 22 of file SignalShapingServiceDUNE35t_service.cc.

  : fInit(false)
{
  // This class is deprecated. See https://cdcvs.fnal.gov/redmine/issues/11777.
  mf::LogInfo("SignalShapingServiceDUNE10kt") << "Deprecated: Consider using SignalShapingServiceDUNE";
  reconfigure(pset);
}

util::SignalShapingServiceDUNE35t::~SignalShapingServiceDUNE35t

(

)

Definition at line 34 of file SignalShapingServiceDUNE35t_service.cc.

{}

Member Function Documentation

template<class T >

void util::SignalShapingServiceDUNE35t::Convolute ( detinfo::DetectorClocksData const &  clockData, unsigned int  channel, std::vector< T > &  func  ) const

inline

Definition at line 179 of file SignalShapingServiceDUNE35t.h.

{
  SignalShaping(channel).Convolute(func);

  //negative number
  int time_offset = FieldResponseTOffset(clockData, channel);
  
  std::vector<T> temp;
  if (time_offset <=0){
    temp.assign(func.begin(),func.begin()-time_offset);
    func.erase(func.begin(),func.begin()-time_offset);
    func.insert(func.end(),temp.begin(),temp.end());
  }else{
    temp.assign(func.end()-time_offset,func.end());
    func.erase(func.end()-time_offset,func.end());
    func.insert(func.begin(),temp.begin(),temp.end());
  }
}

template<class T >

void util::SignalShapingServiceDUNE35t::Deconvolute ( detinfo::DetectorClocksData const &  clockData, unsigned int  channel, std::vector< T > &  func  ) const

inline

Definition at line 202 of file SignalShapingServiceDUNE35t.h.

{
  SignalShaping(channel).Deconvolute(func);

  int time_offset = FieldResponseTOffset(clockData, channel);
  
  std::vector<T> temp;
  if (time_offset <=0){
    temp.assign(func.end()+time_offset,func.end());
    func.erase(func.end()+time_offset,func.end());
    func.insert(func.begin(),temp.begin(),temp.end());
  }else{
    temp.assign(func.begin(),func.begin()+time_offset);
    func.erase(func.begin(),func.begin()+time_offset);
    func.insert(func.end(),temp.begin(),temp.end());    
  }
}

int util::SignalShapingServiceDUNE35t::FieldResponseTOffset	(	detinfo::DetectorClocksData const &	clockData,
		unsigned int const	channel
	)		const

Definition at line 778 of file SignalShapingServiceDUNE35t_service.cc.

{
  art::ServiceHandle<geo::Geometry> geom;
  //geo::SigType_t sigtype = geom->SignalType(channel);
 
  // we need to distinguis between the U and V planes
  geo::View_t view = geom->View(channel);

  double time_offset = 0;

  if(view == geo::kU)
    time_offset = fFieldResponseTOffset.at(0) + fCalibResponseTOffset.at(0);   
  else if(view == geo::kV)
    time_offset = fFieldResponseTOffset.at(1) + fCalibResponseTOffset.at(1);
  else if(view == geo::kZ)
    time_offset = fFieldResponseTOffset.at(2) + fCalibResponseTOffset.at(2); 
  else
    throw cet::exception("SignalShapingServiceDUNE35t")<< "can't determine"
                                                       << " View\n";
 
  return clockData.TPCClock().Ticks(time_offset/1.e3);
  
}

double util::SignalShapingServiceDUNE35t::GetASICGain

(

unsigned int const

channel

)

const

Definition at line 215 of file SignalShapingServiceDUNE35t_service.cc.

{
  art::ServiceHandle<geo::Geometry> geom;
    
  //geo::SigType_t sigtype = geom->SignalType(channel);
  
  // we need to distinguis between the U and V planes
  geo::View_t view = geom->View(channel); 
  
  double gain = 0;
  if(view == geo::kU)
    gain = fASICGainInMVPerFC.at(0);
  else if(view == geo::kV)
    gain = fASICGainInMVPerFC.at(1);
  else if(view == geo::kZ)
    gain = fASICGainInMVPerFC.at(2);
  else
    throw cet::exception("SignalShapingServiceDUNE35t")<< "can't determine"
                                                       << " View\n";
  return gain;
}

double util::SignalShapingServiceDUNE35t::GetDeconNoise

(

unsigned int const

channel

)

const

Definition at line 301 of file SignalShapingServiceDUNE35t_service.cc.

{
  unsigned int plane;
  art::ServiceHandle<geo::Geometry> geom;
  //geo::SigType_t sigtype = geom->SignalType(channel);
  
  // we need to distinguis between the U and V planes
  geo::View_t view = geom->View(channel);
  
  if(view == geo::kU)
    plane = 0;
  else if(view == geo::kV)
    plane = 1;
  else if(view == geo::kZ)
    plane = 2;
  else
    throw cet::exception("SignalShapingServiceDUNE35t")<< "can't determine"
                                                          << " View\n";

  double shapingtime = fShapeTimeConst.at(plane);
  int temp;
  if (shapingtime == 0.5){
    temp = 0;
  }else if (shapingtime == 1.0){
    temp = 1;
  }else if (shapingtime == 2.0){
    temp = 2;
  }else{
    temp = 3;
  }
  auto tempNoise = fNoiseFactVec.at(plane);
  double deconNoise = tempNoise.at(temp);

  deconNoise = deconNoise /4096.*2000./4.7 *6.241*1000/fDeconNorm;
  return deconNoise;
}

double util::SignalShapingServiceDUNE35t::GetDeconNorm ( )

inline

Definition at line 85 of file SignalShapingServiceDUNE35t.h.

{return fDeconNorm;};

std::vector<DoubleVec> util::SignalShapingServiceDUNE35t::GetNoiseFactVec ( )

inline

Definition at line 61 of file SignalShapingServiceDUNE35t.h.

{ return fNoiseFactVec; }

double util::SignalShapingServiceDUNE35t::GetRawNoise

(

unsigned int const

channel

)

const

Definition at line 261 of file SignalShapingServiceDUNE35t_service.cc.

{
  unsigned int plane;
  art::ServiceHandle<geo::Geometry> geom;
  //geo::SigType_t sigtype = geom->SignalType(channel);

  // we need to distinguis between the U and V planes
  geo::View_t view = geom->View(channel);

  if(view == geo::kU)
    plane = 0;
  else if(view == geo::kV)
    plane = 1;
  else if(view == geo::kZ)
    plane = 2;
  else
    throw cet::exception("SignalShapingServiceDUNE35t")<< "can't determine"
                                                          << " View\n";

  double shapingtime = fShapeTimeConst.at(plane);
  double gain = fASICGainInMVPerFC.at(plane);
  int temp;
  if (shapingtime == 0.5){
    temp = 0;
  }else if (shapingtime == 1.0){
    temp = 1;
  }else if (shapingtime == 2.0){
    temp = 2;
  }else{
    temp = 3;
  }
  double rawNoise;

  auto tempNoise = fNoiseFactVec.at(plane);
  rawNoise = tempNoise.at(temp);

  rawNoise *= gain/4.7;
  return rawNoise;
}

double util::SignalShapingServiceDUNE35t::GetShapingTime

(

unsigned int const

channel

)

const

Definition at line 239 of file SignalShapingServiceDUNE35t_service.cc.

{
  art::ServiceHandle<geo::Geometry> geom;
  //geo::SigType_t sigtype = geom->SignalType(channel);

  // we need to distinguis between the U and V planes
  geo::View_t view = geom->View(channel); 

  double shaping_time = 0;

  if(view == geo::kU)
    shaping_time = fShapeTimeConst.at(0);
  else if(view == geo::kV)
     shaping_time = fShapeTimeConst.at(1);
  else if(view == geo::kZ)
    shaping_time = fShapeTimeConst.at(2);
  else
    throw cet::exception("SignalShapingServiceDUNE35t")<< "can't determine"
                                                       << " View\n";
  return shaping_time;
}

void util::SignalShapingServiceDUNE35t::init ( ) const

inlineprivate

Definition at line 93 of file SignalShapingServiceDUNE35t.h.

{const_cast<SignalShapingServiceDUNE35t*>(this)->init();}

void util::SignalShapingServiceDUNE35t::init ( )

private

Definition at line 342 of file SignalShapingServiceDUNE35t_service.cc.

{
  if(!fInit) {
    fInit = true;

    // Do microboone-specific configuration of SignalShaping by providing
    // microboone response and filter functions.

    // Calculate field and electronics response functions.

    // Lazy evaluation in a multi-threaded environment is problematic.
    // The best we can do here is to use the global data available
    // through detector-properties service.
    auto const clockData = art::ServiceHandle<detinfo::DetectorClocksService const>()->DataForJob();
    auto const detProp = art::ServiceHandle<detinfo::DetectorPropertiesService const>()->DataForJob(clockData);
    SetFieldResponse(clockData, detProp);
    SetElectResponse(fShapeTimeConst.at(2),fASICGainInMVPerFC.at(2));

    // Configure convolution kernels.

    fColSignalShaping.AddResponseFunction(fColFieldResponse);
    fColSignalShaping.AddResponseFunction(fElectResponse);
    fColSignalShaping.save_response();
    fColSignalShaping.set_normflag(false);
    //fColSignalShaping.SetPeakResponseTime(0.);

    SetElectResponse(fShapeTimeConst.at(0),fASICGainInMVPerFC.at(0));

    fIndUSignalShaping.AddResponseFunction(fIndUFieldResponse);
    fIndUSignalShaping.AddResponseFunction(fElectResponse);
    fIndUSignalShaping.save_response();
    fIndUSignalShaping.set_normflag(false);
    //fIndSignalShaping.SetPeakResponseTime(0.);

    SetElectResponse(fShapeTimeConst.at(1),fASICGainInMVPerFC.at(1));

    fIndVSignalShaping.AddResponseFunction(fIndVFieldResponse);
    fIndVSignalShaping.AddResponseFunction(fElectResponse);
    fIndVSignalShaping.save_response();
    fIndVSignalShaping.set_normflag(false);
    //fIndSignalShaping.SetPeakResponseTime(0.);
    

    SetResponseSampling(clockData);

    // Calculate filter functions.

    SetFilters(clockData);

    // Configure deconvolution kernels.

    fColSignalShaping.AddFilterFunction(fColFilter);
    fColSignalShaping.CalculateDeconvKernel();

    fIndUSignalShaping.AddFilterFunction(fIndUFilter);
    fIndUSignalShaping.CalculateDeconvKernel();

    fIndVSignalShaping.AddFilterFunction(fIndVFilter);
    fIndVSignalShaping.CalculateDeconvKernel();
  }
}

void util::SignalShapingServiceDUNE35t::reconfigure

(

const fhicl::ParameterSet &

pset

)

Definition at line 40 of file SignalShapingServiceDUNE35t_service.cc.

{
  // Reset initialization flag.

  fInit = false;

  // Reset kernels.

  fColSignalShaping.Reset();
  fIndUSignalShaping.Reset();
  fIndVSignalShaping.Reset();

  // Fetch fcl parameters.

  fNFieldBins = pset.get<int>("FieldBins");
  fCol3DCorrection = pset.get<double>("Col3DCorrection");
  fInd3DCorrection = pset.get<double>("Ind3DCorrection");
  fColFieldRespAmp = pset.get<double>("ColFieldRespAmp");
  fIndUFieldRespAmp = pset.get<double>("IndUFieldRespAmp");
  fIndVFieldRespAmp = pset.get<double>("IndVFieldRespAmp");
  
  fDeconNorm = pset.get<double>("DeconNorm");
  fADCPerPCAtLowestASICGain = pset.get<double>("ADCPerPCAtLowestASICGain");
  fASICGainInMVPerFC = pset.get<std::vector<double> >("ASICGainInMVPerFC");
  fShapeTimeConst = pset.get<std::vector<double> >("ShapeTimeConst");
  fNoiseFactVec =  pset.get<std::vector<DoubleVec> >("NoiseFactVec");

  fInputFieldRespSamplingPeriod = pset.get<double>("InputFieldRespSamplingPeriod");
  
  fFieldResponseTOffset = pset.get<std::vector<double> >("FieldResponseTOffset");
  fCalibResponseTOffset = pset.get<std::vector<double> >("CalibResponseTOffset");

  fUseFunctionFieldShape= pset.get<bool>("UseFunctionFieldShape");
  fUseHistogramFieldShape = pset.get<bool>("UseHistogramFieldShape");

  fGetFilterFromHisto= pset.get<bool>("GetFilterFromHisto");

  fDebugFieldShape       = pset.get<bool>("DebugFieldShape",false);  

  // Construct parameterized collection filter function.
 if(!fGetFilterFromHisto)
 {
  mf::LogInfo("SignalShapingServiceDUNE35t") << "Getting Filter from .fcl file" ;
  std::string colFilt = pset.get<std::string>("ColFilter");
  std::vector<double> colFiltParams =
  pset.get<std::vector<double> >("ColFilterParams");
  fColFilterFunc = new TF1("colFilter", colFilt.c_str());
  for(unsigned int i=0; i<colFiltParams.size(); ++i)
    fColFilterFunc->SetParameter(i, colFiltParams[i]);

  // Construct parameterized induction filter function.

  std::string indUFilt = pset.get<std::string>("IndUFilter");
  std::vector<double> indUFiltParams = pset.get<std::vector<double> >("IndUFilterParams");
  fIndUFilterFunc = new TF1("indUFilter", indUFilt.c_str());
  for(unsigned int i=0; i<indUFiltParams.size(); ++i)
    fIndUFilterFunc->SetParameter(i, indUFiltParams[i]);

  std::string indVFilt = pset.get<std::string>("IndVFilter");
  std::vector<double> indVFiltParams = pset.get<std::vector<double> >("IndVFilterParams");
  fIndVFilterFunc = new TF1("indVFilter", indVFilt.c_str());
  for(unsigned int i=0; i<indVFiltParams.size(); ++i)
    fIndVFilterFunc->SetParameter(i, indVFiltParams[i]);

 }
 else
 {  
   std::string histoname = pset.get<std::string>("FilterHistoName");
   mf::LogInfo("SignalShapingServiceDUNE35t") << " using filter from .root file " ;
   int fNPlanes=3;
   
  // constructor decides if initialized value is a path or an environment variable
  std::string fname;
  cet::search_path sp("FW_SEARCH_PATH");
  sp.find_file(pset.get<std::string>("FilterFunctionFname"), fname);
    
  TFile * in=new TFile(fname.c_str(),"READ");
   for(int i=0;i<fNPlanes;i++){
     TH1D * temp=(TH1D *)in->Get(Form(histoname.c_str(),i));
     fFilterHist[i]=new TH1D(Form(histoname.c_str(),i),Form(histoname.c_str(),i),temp->GetNbinsX(),0,temp->GetNbinsX());
     temp->Copy(*fFilterHist[i]); 
    }
   
   in->Close();
   
 }

 /////////////////////////////////////
 if(fUseFunctionFieldShape)
 {
  std::string colField = pset.get<std::string>("ColFieldShape");
  std::vector<double> colFieldParams =
    pset.get<std::vector<double> >("ColFieldParams");
  fColFieldFunc = new TF1("colField", colField.c_str());
  for(unsigned int i=0; i<colFieldParams.size(); ++i)
    fColFieldFunc->SetParameter(i, colFieldParams[i]);

  // Construct parameterized induction filter function.

  std::string indUField = pset.get<std::string>("IndUFieldShape");
  std::vector<double> indUFieldParams = pset.get<std::vector<double> >("IndUFieldParams");
  fIndUFieldFunc = new TF1("indUField", indUField.c_str());
  for(unsigned int i=0; i<indUFieldParams.size(); ++i)
    fIndUFieldFunc->SetParameter(i, indUFieldParams[i]);
   // Warning, last parameter needs to be multiplied by the FFTSize, in current version of the code,
   
  std::string indVField = pset.get<std::string>("IndVFieldShape");
  std::vector<double> indVFieldParams = pset.get<std::vector<double> >("IndVFieldParams");
  fIndVFieldFunc = new TF1("indVField", indVField.c_str());
  for(unsigned int i=0; i<indVFieldParams.size(); ++i)
    fIndVFieldFunc->SetParameter(i, indVFieldParams[i]);
   // Warning, last parameter needs to be multiplied by the FFTSize, in current version of the code,
   
  } else if ( fUseHistogramFieldShape ) {
    mf::LogInfo("SignalShapingServiceDUNE35t") << " using the field response provided from a .root file " ;
    int fNPlanes = 3;
    
    // constructor decides if initialized value is a path or an environment variable
    std::string fname;   
    cet::search_path sp("FW_SEARCH_PATH");
    sp.find_file( pset.get<std::string>("FieldResponseFname"), fname );
    std::string histoname = pset.get<std::string>("FieldResponseHistoName");

    std::unique_ptr<TFile> fin(new TFile(fname.c_str(), "READ"));
    if ( !fin->IsOpen() ) throw art::Exception( art::errors::NotFound ) << "Could not find the field response file " << fname << "!" << std::endl;

    std::string iPlane[3] = { "U", "V", "Y" };

    for ( int i = 0; i < fNPlanes; i++ ) {
      TString iHistoName = Form( "%s_%s", histoname.c_str(), iPlane[i].c_str());
      TH1F *temp = (TH1F*) fin->Get( iHistoName );  
      if ( !temp ) throw art::Exception( art::errors::NotFound ) << "Could not find the field response histogram " << iHistoName << std::endl;
      if ( temp->GetNbinsX() > fNFieldBins ) throw art::Exception( art::errors::InvalidNumber ) << "FieldBins should always be larger than or equal to the number of the bins in the input histogram!" << std::endl;
      
      fFieldResponseHist[i] = new TH1F( iHistoName, iHistoName, temp->GetNbinsX(), temp->GetBinLowEdge(1), temp->GetBinLowEdge( temp->GetNbinsX() + 1) );
      temp->Copy(*fFieldResponseHist[i]);
    }
    
    fin->Close();
  }
}

void util::SignalShapingServiceDUNE35t::SetElectResponse ( double  shapingtime, double  gain  )

private

Definition at line 528 of file SignalShapingServiceDUNE35t_service.cc.

{
  // Get services.

  art::ServiceHandle<geo::Geometry> geo;
  art::ServiceHandle<util::LArFFT> fft;

  MF_LOG_DEBUG("SignalShapingDUNE35t") << "Setting DUNE35t electronics response function...";

  int nticks = fft->FFTSize();
  fElectResponse.resize(nticks, 0.);
  std::vector<double> time(nticks,0.);

  //Gain and shaping time variables from fcl file:    
  double Ao = 1.0;
  double To = shapingtime;  //peaking time
    
  // this is actually sampling time, in ns
  // mf::LogInfo("SignalShapingDUNE35t") << "Check sampling intervals: " 
  //                                  << fSampleRate << " ns" 
  //                                  << "Check number of samples: " << fNTicks;

  // The following sets the microboone electronics response function in 
  // time-space. Function comes from BNL SPICE simulation of DUNE35t 
  // electronics. SPICE gives the electronics transfer function in 
  // frequency-space. The inverse laplace transform of that function 
  // (in time-space) was calculated in Mathematica and is what is being 
  // used below. Parameters Ao and To are cumulative gain/timing parameters 
  // from the full (ASIC->Intermediate amp->Receiver->ADC) electronics chain. 
  // They have been adjusted to make the SPICE simulation to match the 
  // actual electronics response. Default params are Ao=1.4, To=0.5us. 
  double max = 0;
  
  for(size_t i = 0; i < fElectResponse.size(); ++i){

    //convert time to microseconds, to match fElectResponse[i] definition
    time[i] = (1.*i)*fInputFieldRespSamplingPeriod *1e-3; 
    fElectResponse[i] = 
      4.31054*exp(-2.94809*time[i]/To)*Ao - 2.6202*exp(-2.82833*time[i]/To)*cos(1.19361*time[i]/To)*Ao
      -2.6202*exp(-2.82833*time[i]/To)*cos(1.19361*time[i]/To)*cos(2.38722*time[i]/To)*Ao
      +0.464924*exp(-2.40318*time[i]/To)*cos(2.5928*time[i]/To)*Ao
      +0.464924*exp(-2.40318*time[i]/To)*cos(2.5928*time[i]/To)*cos(5.18561*time[i]/To)*Ao
      +0.762456*exp(-2.82833*time[i]/To)*sin(1.19361*time[i]/To)*Ao
      -0.762456*exp(-2.82833*time[i]/To)*cos(2.38722*time[i]/To)*sin(1.19361*time[i]/To)*Ao
      +0.762456*exp(-2.82833*time[i]/To)*cos(1.19361*time[i]/To)*sin(2.38722*time[i]/To)*Ao
      -2.6202*exp(-2.82833*time[i]/To)*sin(1.19361*time[i]/To)*sin(2.38722*time[i]/To)*Ao 
      -0.327684*exp(-2.40318*time[i]/To)*sin(2.5928*time[i]/To)*Ao + 
      +0.327684*exp(-2.40318*time[i]/To)*cos(5.18561*time[i]/To)*sin(2.5928*time[i]/To)*Ao
      -0.327684*exp(-2.40318*time[i]/To)*cos(2.5928*time[i]/To)*sin(5.18561*time[i]/To)*Ao
      +0.464924*exp(-2.40318*time[i]/To)*sin(2.5928*time[i]/To)*sin(5.18561*time[i]/To)*Ao;

    if(fElectResponse[i] > max) max = fElectResponse[i];
    
  }// end loop over time buckets
    

  MF_LOG_DEBUG("SignalShapingDUNE35t") << " Done.";

 //normalize fElectResponse[i], before the convolution   
  
   for(auto& element : fElectResponse){
    element /= max;
    element *= fADCPerPCAtLowestASICGain * 1.60217657e-7;
    element *= gain / 4.7;
   }
  
  return;

}

void util::SignalShapingServiceDUNE35t::SetFieldResponse ( detinfo::DetectorClocksData const &  clockData, detinfo::DetectorPropertiesData const &  detProp  )

private

in cm

number of bins //KP

Definition at line 407 of file SignalShapingServiceDUNE35t_service.cc.

{
  // Get services.

  art::ServiceHandle<geo::Geometry> geo;

  // Get plane pitch.
 
  double xyz1[3] = {0.};
  double xyz2[3] = {0.};
  double xyzl[3] = {0.};
  // should always have at least 2 planes
  geo->Plane(0).LocalToWorld(xyzl, xyz1);
  geo->Plane(1).LocalToWorld(xyzl, xyz2);

  // this assumes all planes are equidistant from each other,
  // probably not a bad assumption
  double pitch = xyz2[0] - xyz1[0]; ///in cm

  fColFieldResponse.resize(fNFieldBins, 0.);
  fIndUFieldResponse.resize(fNFieldBins, 0.);
  fIndVFieldResponse.resize(fNFieldBins, 0.);

  // set the response for the collection plane first
  // the first entry is 0

  double driftvelocity=detProp.DriftVelocity()/1000.;
  int nbinc = TMath::Nint(fCol3DCorrection*(fabs(pitch))/(driftvelocity*sampling_rate(clockData))); ///number of bins //KP
  double integral = 0.;  
  ////////////////////////////////////////////////////
   if(fUseFunctionFieldShape)
  {
  art::ServiceHandle<util::LArFFT> fft;
  int signalSize = fft->FFTSize();
  std::vector<double> ramp(signalSize);
   // TComplex kernBin;
   // int size = signalSize/2;
   // int bin=0;
    //std::vector<TComplex> freqSig(size+1);
  std::vector<double> bipolar(signalSize);
    
    
  fColFieldResponse.resize(signalSize, 0.);
  fIndUFieldResponse.resize(signalSize, 0.);
  fIndVFieldResponse.resize(signalSize, 0.);
   
  // Hardcoding. Bad. Temporary hopefully.
  fIndUFieldFunc->SetParameter(4,fIndUFieldFunc->GetParameter(4)*signalSize);
  fIndVFieldFunc->SetParameter(4,fIndVFieldFunc->GetParameter(4)*signalSize);
  
  //double integral = 0.;
    for(int i = 0; i < signalSize; i++) {
          ramp[i]=fColFieldFunc->Eval(i);
          fColFieldResponse[i]=ramp[i];
          integral += fColFieldResponse[i];
          // rampc->Fill(i,ramp[i]);
          bipolar[i]=fIndUFieldFunc->Eval(i);
          fIndUFieldResponse[i]=bipolar[i];
          bipolar[i]=fIndVFieldFunc->Eval(i);
          fIndVFieldResponse[i]=bipolar[i];
          // bipol->Fill(i,bipolar[i]);
    }
     
   for(int i = 0; i < signalSize; ++i){
      fColFieldResponse[i] *= fColFieldRespAmp/integral;
     }
      
    //this might be not necessary if the function definition is not defined in the middle of the signal range  
    fft->ShiftData(fIndUFieldResponse,signalSize/2.0);
  } else if ( fUseHistogramFieldShape ) {
    
    // Ticks in nanosecond
    // Calculate the normalization of the collection plane
    for ( int ibin = 1; ibin <= fFieldResponseHist[2]->GetNbinsX(); ibin++ )
      integral += fFieldResponseHist[2]->GetBinContent( ibin );   

    // Induction plane
    for ( int ibin = 1; ibin <= fFieldResponseHist[0]->GetNbinsX(); ibin++ )
      fIndUFieldResponse[ibin-1] = fIndUFieldRespAmp*fFieldResponseHist[0]->GetBinContent( ibin )/integral;

    for ( int ibin = 1; ibin <= fFieldResponseHist[1]->GetNbinsX(); ibin++ )
      fIndVFieldResponse[ibin-1] = fIndVFieldRespAmp*fFieldResponseHist[1]->GetBinContent( ibin )/integral;

    // Collection plane
    for ( int ibin = 1; ibin <= fFieldResponseHist[2]->GetNbinsX(); ibin++ )
      fColFieldResponse[ibin-1] = fColFieldRespAmp*fFieldResponseHist[2]->GetBinContent( ibin )/integral;
  } else
  {
    //////////////////////////////////////////////////
    mf::LogInfo("SignalShapingServiceDUNE35t") << " using the old field shape " ;
    double integral = 0.;
    for(int i = 1; i < nbinc; ++i){
      fColFieldResponse[i] = fColFieldResponse[i-1] + 1.0;
      integral += fColFieldResponse[i];
    }
    
    for(int i = 0; i < nbinc; ++i){
      fColFieldResponse[i] *= fColFieldRespAmp/integral;
    }
    
    // now the induction plane
    
    int nbini = TMath::Nint(fInd3DCorrection*(fabs(pitch))/(driftvelocity*sampling_rate(clockData)));//KP
    for(int i = 0; i < nbini; ++i){
      fIndUFieldResponse[i] = fIndUFieldRespAmp/(1.*nbini);
      fIndUFieldResponse[nbini+i] = -fIndUFieldRespAmp/(1.*nbini);
    }

    for(int i = 0; i < nbini; ++i){
      fIndVFieldResponse[i] = fIndVFieldRespAmp/(1.*nbini);
      fIndVFieldResponse[nbini+i] = -fIndVFieldRespAmp/(1.*nbini);
    }

  }
  
  return;
}

void util::SignalShapingServiceDUNE35t::SetFilters ( detinfo::DetectorClocksData const &  clockData )

private

Definition at line 601 of file SignalShapingServiceDUNE35t_service.cc.

{
  // Get services.

  art::ServiceHandle<util::LArFFT> fft;

  double ts = sampling_rate(clockData);
  int n = fft->FFTSize() / 2;

  // Calculate collection filter.

  fColFilter.resize(n+1);
  fIndUFilter.resize(n+1);
  fIndVFilter.resize(n+1);
  
  if(!fGetFilterFromHisto)
  {
  fColFilterFunc->SetRange(0, double(n));

  for(int i=0; i<=n; ++i) {
    double freq = 500. * i / (ts * n);      // Cycles / microsecond.
    double f = fColFilterFunc->Eval(freq);
    fColFilter[i] = TComplex(f, 0.);
  }
  

  // Calculate induction filter.

 
  fIndUFilterFunc->SetRange(0, double(n));

  for(int i=0; i<=n; ++i) {
    double freq = 500. * i / (ts * n);      // Cycles / microsecond.
    double f = fIndUFilterFunc->Eval(freq);
    fIndUFilter[i] = TComplex(f, 0.);
    }
  
  fIndVFilterFunc->SetRange(0, double(n));

  for(int i=0; i<=n; ++i) {
    double freq = 500. * i / (ts * n);      // Cycles / microsecond.
    double f = fIndVFilterFunc->Eval(freq);
    fIndVFilter[i] = TComplex(f, 0.);
  } 

  }
  else
  {
    
    for(int i=0; i<=n; ++i) {
      double f = fFilterHist[2]->GetBinContent(i);  // hardcoded plane numbers. Bad. To change later.
      fColFilter[i] = TComplex(f, 0.);
      double g = fFilterHist[1]->GetBinContent(i);
      fIndVFilter[i] = TComplex(g, 0.);
      double h = fFilterHist[0]->GetBinContent(i);
      fIndUFilter[i] = TComplex(h, 0.); 
    }
  }
  
  // fIndUSignalShaping.AddFilterFunction(fIndUFilter);
  // fIndVSignalShaping.AddFilterFunction(fIndVFilter);
  // fColSignalShaping.AddFilterFunction(fColFilter);
  
}

void util::SignalShapingServiceDUNE35t::SetResponseSampling ( detinfo::DetectorClocksData const &  clockData )

private

VELOCITY-OUT ... comment out kDVel usage here

VELOCITY-OUT ... comment out kDVel usage here

VELOCITY-OUT ... comment out kDVel usage here

Definition at line 670 of file SignalShapingServiceDUNE35t_service.cc.

{
  // Get services
  art::ServiceHandle<geo::Geometry> geo;
  art::ServiceHandle<util::LArFFT> fft;

  // Operation permitted only if output of rebinning has a larger bin size
  if( fInputFieldRespSamplingPeriod > sampling_rate(clockData) )
    throw cet::exception(__FUNCTION__) << "\033[93m"
                                       << "Invalid operation: cannot rebin to a more finely binned vector!"
                                       << "\033[00m" << std::endl;

  int nticks = fft->FFTSize();
  std::vector<double> SamplingTime( nticks, 0. );
  for ( int itime = 0; itime < nticks; itime++ ) {
    SamplingTime[itime] = (1.*itime) * sampling_rate(clockData);
    /// VELOCITY-OUT ... comment out kDVel usage here
    //SamplingTime[itime] = (1.*itime) * sampling_rate(clockData) / kDVel;
  }

  // Sampling
  for ( int iplane = 0; iplane < 3; iplane++ ) {
    const std::vector<double>* pResp;
    switch ( iplane ) {
    case 0: pResp = &(fIndUSignalShaping.Response_save()); break;
    case 1: pResp = &(fIndVSignalShaping.Response_save()); break;
    default: pResp = &(fColSignalShaping.Response_save()); break;
    }

    std::vector<double> SamplingResp(nticks , 0. );
    
    
    int nticks_input = pResp->size();
    std::vector<double> InputTime(nticks_input, 0. );
    for ( int itime = 0; itime < nticks_input; itime++ ) {
      InputTime[itime] = (1.*itime) * fInputFieldRespSamplingPeriod;
    }
    
   
    /*
      Much more sophisticated approach using a linear (trapezoidal) interpolation 
      Current default!
    */
    int SamplingCount = 0;    
    for ( int itime = 0; itime < nticks; itime++ ) {
      int low = -1, up = -1;
      for ( int jtime = 0; jtime < nticks_input; jtime++ ) {
        if ( InputTime[jtime] == SamplingTime[itime] ) {
          SamplingResp[itime] = (*pResp)[jtime];
          /// VELOCITY-OUT ... comment out kDVel usage here
          //SamplingResp[itime] = kDVel * (*pResp)[jtime];
          SamplingCount++;
          break;
        } else if ( InputTime[jtime] > SamplingTime[itime] ) {
          low = jtime - 1;
          up = jtime;
          SamplingResp[itime] = (*pResp)[low] + ( SamplingTime[itime] - InputTime[low] ) * ( (*pResp)[up] - (*pResp)[low] ) / ( InputTime[up] - InputTime[low] );
          /// VELOCITY-OUT ... comment out kDVel usage here
          //SamplingResp[itime] *= kDVel;
          SamplingCount++;
          break;
        } else {
          SamplingResp[itime] = 0.;
        }
      } // for ( int jtime = 0; jtime < nticks; jtime++ )
    } // for ( int itime = 0; itime < nticks; itime++ )
    SamplingResp.resize( SamplingCount, 0.);    

  

  
    switch ( iplane ) {
    case 0: fIndUSignalShaping.AddResponseFunction( SamplingResp, true ); break;
    case 1: fIndVSignalShaping.AddResponseFunction( SamplingResp, true ); break;
    default: fColSignalShaping.AddResponseFunction( SamplingResp, true ); break;
    }

    if (fDebugFieldShape){
      if (iplane == 0){
        std::ofstream outfile("resutest.txt");
        for (size_t i = 0; i<SamplingResp.size(); ++i){
          outfile<<i<<" "<<SamplingResp[i]<<std::endl;
        }
        outfile.close();
      }
      else if (iplane==1){
        std::ofstream outfile("resvtest.txt");
        for (size_t i = 0; i<SamplingResp.size(); ++i){
          outfile<<i<<" "<<SamplingResp[i]<<std::endl;
        }
        outfile.close();
      }
      else if (iplane==2){
        std::ofstream outfile("resztest.txt");
        for (size_t i = 0; i<SamplingResp.size(); ++i){
          outfile<<i<<" "<<SamplingResp[i]<<std::endl;
        }
        outfile.close();
      }
    }   

  } // for ( int iplane = 0; iplane < fNPlanes; iplane++ )

  return;
}

const util::SignalShaping & util::SignalShapingServiceDUNE35t::SignalShaping

(

unsigned int

channel

)

const

Definition at line 186 of file SignalShapingServiceDUNE35t_service.cc.

{
  if(!fInit)
    init();

  // Figure out plane type.

  art::ServiceHandle<geo::Geometry> geom;
  //geo::SigType_t sigtype = geom->SignalType(channel);

  // we need to distinguis between the U and V planes
  geo::View_t view = geom->View(channel); 

  // Return appropriate shaper.

  if(view == geo::kU)
    return fIndUSignalShaping;
  else if (view == geo::kV)
    return fIndVSignalShaping;
  else if(view == geo::kZ)
    return fColSignalShaping;
  else
    throw cet::exception("SignalShapingServiceDUNE35t")<< "can't determine"
                                                          << " View\n";
                                                          
return fColSignalShaping;
}

Member Data Documentation

double util::SignalShapingServiceDUNE35t::fADCPerPCAtLowestASICGain

private

Pulse amplitude gain for a 1 pc charge impulse after convoluting it the with field and electronics response with the lowest ASIC gain setting of 4.7 mV/fC.

Definition at line 132 of file SignalShapingServiceDUNE35t.h.

std::vector<double> util::SignalShapingServiceDUNE35t::fASICGainInMVPerFC

private

Definition at line 135 of file SignalShapingServiceDUNE35t.h.

std::vector<double> util::SignalShapingServiceDUNE35t::fCalibResponseTOffset

private

Definition at line 128 of file SignalShapingServiceDUNE35t.h.

double util::SignalShapingServiceDUNE35t::fCol3DCorrection

private

correction factor to account for 3D path of electrons thru wires

Definition at line 119 of file SignalShapingServiceDUNE35t.h.

TF1* util::SignalShapingServiceDUNE35t::fColFieldFunc

private

Parameterized collection field shape function.

Definition at line 145 of file SignalShapingServiceDUNE35t.h.

double util::SignalShapingServiceDUNE35t::fColFieldRespAmp

private

amplitude of response to field

Definition at line 123 of file SignalShapingServiceDUNE35t.h.

std::vector<double> util::SignalShapingServiceDUNE35t::fColFieldResponse

private

Definition at line 159 of file SignalShapingServiceDUNE35t.h.

std::vector<TComplex> util::SignalShapingServiceDUNE35t::fColFilter

private

Definition at line 169 of file SignalShapingServiceDUNE35t.h.

TF1* util::SignalShapingServiceDUNE35t::fColFilterFunc

private

Parameterized collection filter function.

Definition at line 138 of file SignalShapingServiceDUNE35t.h.

util::SignalShaping util::SignalShapingServiceDUNE35t::fColSignalShaping

private

Definition at line 153 of file SignalShapingServiceDUNE35t.h.

bool util::SignalShapingServiceDUNE35t::fDebugFieldShape

private

Save field response shape to a text file.

Definition at line 173 of file SignalShapingServiceDUNE35t.h.

double util::SignalShapingServiceDUNE35t::fDeconNorm

private

Definition at line 131 of file SignalShapingServiceDUNE35t.h.

std::vector<double> util::SignalShapingServiceDUNE35t::fElectResponse

private

Definition at line 165 of file SignalShapingServiceDUNE35t.h.

TH1F* util::SignalShapingServiceDUNE35t::fFieldResponseHist[3]

private

Histogram used to hold the field response, hardcoded for the time being.

Definition at line 148 of file SignalShapingServiceDUNE35t.h.

std::vector<double> util::SignalShapingServiceDUNE35t::fFieldResponseTOffset

private

Time offset for field response in ns.

Definition at line 127 of file SignalShapingServiceDUNE35t.h.

TH1D* util::SignalShapingServiceDUNE35t::fFilterHist[3]

private

Histogram used to hold the collection filter, hardcoded for the time being.

Definition at line 149 of file SignalShapingServiceDUNE35t.h.

bool util::SignalShapingServiceDUNE35t::fGetFilterFromHisto

private

Flag that allows to use a filter function from a histogram instead of the functional dependency.

Definition at line 144 of file SignalShapingServiceDUNE35t.h.

double util::SignalShapingServiceDUNE35t::fInd3DCorrection

private

correction factor to account for 3D path of electrons thru wires

Definition at line 121 of file SignalShapingServiceDUNE35t.h.

TF1* util::SignalShapingServiceDUNE35t::fIndUFieldFunc

private

Parameterized induction field shape function.

Definition at line 146 of file SignalShapingServiceDUNE35t.h.

double util::SignalShapingServiceDUNE35t::fIndUFieldRespAmp

private

amplitude of response to field

Definition at line 124 of file SignalShapingServiceDUNE35t.h.

std::vector<double> util::SignalShapingServiceDUNE35t::fIndUFieldResponse

private

Definition at line 160 of file SignalShapingServiceDUNE35t.h.

std::vector<TComplex> util::SignalShapingServiceDUNE35t::fIndUFilter

private

Definition at line 170 of file SignalShapingServiceDUNE35t.h.

TF1* util::SignalShapingServiceDUNE35t::fIndUFilterFunc

private

Parameterized induction filter function.

Definition at line 139 of file SignalShapingServiceDUNE35t.h.

util::SignalShaping util::SignalShapingServiceDUNE35t::fIndUSignalShaping

private

Definition at line 154 of file SignalShapingServiceDUNE35t.h.

TF1* util::SignalShapingServiceDUNE35t::fIndVFieldFunc

private

Parameterized induction field shape function.

Definition at line 147 of file SignalShapingServiceDUNE35t.h.

double util::SignalShapingServiceDUNE35t::fIndVFieldRespAmp

private

amplitude of response to field

Definition at line 125 of file SignalShapingServiceDUNE35t.h.

std::vector<double> util::SignalShapingServiceDUNE35t::fIndVFieldResponse

private

Definition at line 161 of file SignalShapingServiceDUNE35t.h.

std::vector<TComplex> util::SignalShapingServiceDUNE35t::fIndVFilter

private

Definition at line 171 of file SignalShapingServiceDUNE35t.h.

TF1* util::SignalShapingServiceDUNE35t::fIndVFilterFunc

private

Parameterized induction filter function.

Definition at line 140 of file SignalShapingServiceDUNE35t.h.

util::SignalShaping util::SignalShapingServiceDUNE35t::fIndVSignalShaping

private

Definition at line 155 of file SignalShapingServiceDUNE35t.h.

bool util::SignalShapingServiceDUNE35t::fInit

private

Initialization flag.

Definition at line 109 of file SignalShapingServiceDUNE35t.h.

double util::SignalShapingServiceDUNE35t::fInputFieldRespSamplingPeriod

private

Sampling period in the input field response.

Definition at line 129 of file SignalShapingServiceDUNE35t.h.

int util::SignalShapingServiceDUNE35t::fNFieldBins

private

number of bins for field response

Definition at line 118 of file SignalShapingServiceDUNE35t.h.

std::vector<DoubleVec> util::SignalShapingServiceDUNE35t::fNoiseFactVec

private

Definition at line 133 of file SignalShapingServiceDUNE35t.h.

std::vector<double> util::SignalShapingServiceDUNE35t::fShapeTimeConst

private

time constants for exponential shaping

Definition at line 137 of file SignalShapingServiceDUNE35t.h.

bool util::SignalShapingServiceDUNE35t::fUseFunctionFieldShape

private

Flag that allows to use a parameterized field response instead of the hardcoded version.

Definition at line 142 of file SignalShapingServiceDUNE35t.h.

bool util::SignalShapingServiceDUNE35t::fUseHistogramFieldShape

private

Flag that turns on field response shapes from histograms.

Definition at line 143 of file SignalShapingServiceDUNE35t.h.

---

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

• dunecore/dunecore/Utilities/SignalShapingServiceDUNE35t.h
• dunecore/dunecore/Utilities/SignalShapingServiceDUNE35t_service.cc