CalibrationAlgorithm Class Reference

Algorithm for calibration runs. More...

#include <CalibrationAlgorithm.h>

Inheritance diagram for CalibrationAlgorithm:

Public Member Functions
	CalibrationAlgorithm (const edm::ParameterSet &pset, CalibrationAnalysis *const )
const Histo &	histo (int i) const
	~CalibrationAlgorithm () override
Public Member Functions inherited from CommissioningAlgorithm
void	analysis (const std::vector< TH1 * > &)
	CommissioningAlgorithm (CommissioningAnalysis *const )
	CommissioningAlgorithm ()
virtual	~CommissioningAlgorithm ()

Private Member Functions
void	analyse () override
	CalibrationAlgorithm ()
void	correctDistribution (TH1 *) const
void	extract (const std::vector< TH1 * > &) override
TF1 *	fitPulse (TH1 *, float rangeLow=0, float rangeHigh=-1)
float	maximum (TH1 *)
float	turnOn (TH1 *)

Private Attributes
CalibrationAnalysis *	cal_
TF1 *	deconv_fitter_
Histo	histo_ [32]
TF1 *	peak_fitter_

Additional Inherited Members
Public Types inherited from CommissioningAlgorithm
typedef std::pair< TH1 *, std::string >	Histo
Protected Member Functions inherited from CommissioningAlgorithm
CommissioningAnalysis *const	anal () const
uint32_t	extractFedKey (const TH1 *const )

Detailed Description

Algorithm for calibration runs.

Author

C. Delaere

Definition at line 18 of file CalibrationAlgorithm.h.

Constructor & Destructor Documentation

CalibrationAlgorithm::CalibrationAlgorithm	(	const edm::ParameterSet &	pset,
		CalibrationAnalysis * const	anal
	)

Definition at line 20 of file CalibrationAlgorithm.cc.

References deconv_fitter_, fdeconv_convoluted(), fpeak_convoluted(), and peak_fitter_.

  : CommissioningAlgorithm(anal),
    deconv_fitter_(nullptr),
    peak_fitter_(nullptr),
    cal_(nullptr)
{
  deconv_fitter_ = new TF1("deconv_fitter",fdeconv_convoluted,-50,50,5);
  deconv_fitter_->FixParameter(0,0);
  deconv_fitter_->SetParLimits(1,-100,0);
  deconv_fitter_->SetParLimits(2,0,200);
  deconv_fitter_->SetParLimits(3,5,100);
  deconv_fitter_->FixParameter(3,50);
  deconv_fitter_->SetParLimits(4,0,50);
  deconv_fitter_->SetParameters(0.,-10,0.96,50,20);
  peak_fitter_ = new TF1("peak_fitter",fpeak_convoluted,-50,50,5);
  peak_fitter_->FixParameter(0,0);
  peak_fitter_->SetParLimits(1,-100,0);
  peak_fitter_->SetParLimits(2,0,400);
  peak_fitter_->SetParLimits(3,5,100);
  peak_fitter_->FixParameter(3,50);
  peak_fitter_->SetParLimits(4,0,50);
  peak_fitter_->SetParameters(0.,-10,0.96,50,20);
}

CalibrationAlgorithm::~CalibrationAlgorithm ( )

inlineoverride

Definition at line 24 of file CalibrationAlgorithm.h.

{;}

CalibrationAlgorithm::CalibrationAlgorithm ( )

inlineprivate

Definition at line 30 of file CalibrationAlgorithm.h.

References analyse(), correctDistribution(), extract(), fitPulse(), maximum(), and turnOn().

{;}

Member Function Documentation

void CalibrationAlgorithm::analyse ( )

overrideprivatevirtual

Performs histogram anaylsis.

Implements CommissioningAlgorithm.

Definition at line 101 of file CalibrationAlgorithm.cc.

References CalibrationAnalysis::amplitude_, cal_, CalibrationAnalysis::calchan_, CalibrationAnalysis::chi2_, correctDistribution(), plotBeamSpotDB::first, fitPulse(), histo_, mps_fire::i, CalibrationAnalysis::isScan_, CalibrationAnalysis::max_amplitude_, CalibrationAnalysis::max_chi2_, CalibrationAnalysis::max_riseTime_, CalibrationAnalysis::max_smearing_, CalibrationAnalysis::max_tail_, CalibrationAnalysis::max_timeConstant_, maximum(), CalibrationAnalysis::mean_amplitude_, CalibrationAnalysis::mean_chi2_, CalibrationAnalysis::mean_riseTime_, CalibrationAnalysis::mean_smearing_, CalibrationAnalysis::mean_tail_, CalibrationAnalysis::mean_timeConstant_, CalibrationAnalysis::min_amplitude_, CalibrationAnalysis::min_chi2_, CalibrationAnalysis::min_riseTime_, CalibrationAnalysis::min_smearing_, CalibrationAnalysis::min_tail_, CalibrationAnalysis::min_timeConstant_, sistrip::mlCommissioning_, CalibrationAnalysis::riseTime_, Gflash::Rmax, Gflash::Rmin, CalibrationAnalysis::smearing_, CalibrationAnalysis::spread_amplitude_, CalibrationAnalysis::spread_chi2_, CalibrationAnalysis::spread_riseTime_, CalibrationAnalysis::spread_smearing_, CalibrationAnalysis::spread_tail_, CalibrationAnalysis::spread_timeConstant_, mathSSE::sqrt(), AlCaHLTBitMon_QueryRunRegistry::string, digi_MixPreMix_cfi::strip, CalibrationAnalysis::tail_, CalibrationAnalysis::timeConstant_, fftjetcommon_cfi::title, OfflinePixel3DPrimaryVertices_cfi::Tmin, and turnOn().

Referenced by CalibrationAlgorithm().

                                   { 

  if ( !cal_ ) {
    edm::LogWarning(mlCommissioning_)
      << "[CalibrationAlgorithm::" << __func__ << "]"
      << " NULL pointer to derived Analysis object!";
    return; 
  }

  float Amean[2]   = {0.,0.};
  float Amin[2]    = {2000.,2000.};
  float Amax[2]    = {0.,0.};
  float Aspread[2] = {0.,0.};
  float Tmean[2]   = {0.,0.};
  float Tmin[2]    = {2000.,2000.};
  float Tmax[2]    = {0.,0.};
  float Tspread[2] = {0.,0.};
  float Rmean[2]   = {0.,0.};
  float Rmin[2]    = {2000.,2000.};
  float Rmax[2]    = {0.,0.};
  float Rspread[2] = {0.,0.};
  float Cmean[2]   = {0.,0.};
  float Cmin[2]    = {2000.,2000.};
  float Cmax[2]    = {0.,0.};
  float Cspread[2] = {0.,0.};
  float Smean[2]   = {0.,0.};
  float Smin[2]    = {2000.,2000.};
  float Smax[2]    = {0.,0.};
  float Sspread[2] = {0.,0.};
  float Kmean[2]   = {0.,0.};
  float Kmin[2]    = {2000000.,2000000.};
  float Kmax[2]    = {0.,0.};
  float Kspread[2] = {0.,0.};
 
  unsigned int upperLimit = cal_->isScan_ ? 2 : 32;
  float nStrips = cal_->isScan_ ? 1. : 16.;
  for(unsigned int i=0;i<upperLimit;++i) {
    if ( !histo_[i].first ) {
      edm::LogWarning(mlCommissioning_) 
    << " NULL pointer to histogram " << i << "!";
      return;
    }
     
    // determine which APV and strip is being looked at.
    std::string title = histo_[i].first->GetName();
    int apv = 0;
    int strip = 0;
    if(title.find("STRIP")!=std::string::npos && title.find("Apv")!=std::string::npos) {
      strip = atoi(title.c_str()+title.find("STRIP")+6);
      apv = (atoi(title.c_str()+title.find("Apv")+3))%2;
    } else {
      strip = (atoi(title.c_str()+title.find_last_of("_")+1))%16;
      apv = (atoi(title.c_str()+title.find_last_of("_")+1))/16;
      if(title.find("Apv")!=std::string::npos) {
        apv = (atoi(title.c_str()+title.find("Apv")+3))%2;
    strip = strip*8 + cal_->calchan_;
      } else {
        edm::LogWarning(mlCommissioning_) 
      << " Malformed histogram title! Strip/APV not retreived: " 
      << title;
      }
    }
    
    // rescale the plot
    correctDistribution(histo_[i].first);
    
    // amplitude
    cal_->amplitude_[apv][strip] = histo_[i].first->GetMaximum();
    
    // rise time
    cal_->riseTime_[apv][strip] = maximum(histo_[i].first) - turnOn(histo_[i].first);
    
    // tail 125 ns after the maximum
    int lastBin = histo_[i].first->FindBin(histo_[i].first->GetBinCenter(histo_[i].first->GetMaximumBin())+125);
    if(lastBin>histo_[i].first->GetNbinsX()-4) lastBin = histo_[i].first->GetNbinsX()-4;
     if(histo_[i].first->GetMaximum()!=0)
       cal_->tail_[apv][strip] = 100*histo_[i].first->GetBinContent(lastBin)/histo_[i].first->GetMaximum();
     else
       cal_->tail_[apv][strip] = 100;
  
    // perform the fit for the next quantities
    TF1* fit = fitPulse(histo_[i].first);
    
    // time constant
    cal_->timeConstant_[apv][strip] = fit->GetParameter(3);
  
    // smearing
    cal_->smearing_[apv][strip] = fit->GetParameter(4);
  
    // chi2
    cal_->chi2_[apv][strip] = fit->GetChisquare();
    
    //compute mean, max, min, spread
    Amean[apv] += cal_->amplitude_[apv][strip]/nStrips;
    Amin[apv] = Amin[apv]<cal_->amplitude_[apv][strip] ? Amin[apv] : cal_->amplitude_[apv][strip];
    Amax[apv] = Amax[apv]>cal_->amplitude_[apv][strip] ? Amax[apv] : cal_->amplitude_[apv][strip];
    Aspread[apv] += cal_->amplitude_[apv][strip]*cal_->amplitude_[apv][strip]/nStrips;
    Tmean[apv] += cal_->tail_[apv][strip]/nStrips;
    Tmin[apv] = Tmin[apv]<cal_->tail_[apv][strip] ? Tmin[apv] : cal_->tail_[apv][strip];
    Tmax[apv] = Tmax[apv]>cal_->tail_[apv][strip] ? Tmax[apv] : cal_->tail_[apv][strip];
    Tspread[apv] += cal_->tail_[apv][strip]*cal_->tail_[apv][strip]/nStrips;
    Rmean[apv] += cal_->riseTime_[apv][strip]/nStrips;
    Rmin[apv] = Rmin[apv]<cal_->riseTime_[apv][strip] ? Rmin[apv] : cal_->riseTime_[apv][strip];
    Rmax[apv] = Rmax[apv]>cal_->riseTime_[apv][strip] ? Rmax[apv] : cal_->riseTime_[apv][strip];
    Rspread[apv] += cal_->riseTime_[apv][strip]*cal_->riseTime_[apv][strip]/nStrips;
    Cmean[apv] += cal_->timeConstant_[apv][strip]/nStrips;
    Cmin[apv] = Cmin[apv]<cal_->timeConstant_[apv][strip] ? Cmin[apv] : cal_->timeConstant_[apv][strip];
    Cmax[apv] = Cmax[apv]>cal_->timeConstant_[apv][strip] ? Cmax[apv] : cal_->timeConstant_[apv][strip];
    Cspread[apv] += cal_->timeConstant_[apv][strip]*cal_->timeConstant_[apv][strip]/nStrips;
    Smean[apv] += cal_->smearing_[apv][strip]/nStrips;
    Smin[apv] = Smin[apv]<cal_->smearing_[apv][strip] ? Smin[apv] : cal_->smearing_[apv][strip];
    Smax[apv] = Smax[apv]>cal_->smearing_[apv][strip] ? Smax[apv] : cal_->smearing_[apv][strip];
    Sspread[apv] += cal_->smearing_[apv][strip]*cal_->smearing_[apv][strip]/nStrips;
    Kmean[apv] += cal_->chi2_[apv][strip]/nStrips;
    Kmin[apv] = Kmin[apv]<cal_->chi2_[apv][strip] ? Kmin[apv] : cal_->chi2_[apv][strip];
    Kmax[apv] = Kmax[apv]>cal_->chi2_[apv][strip] ? Kmax[apv] : cal_->chi2_[apv][strip];
    Kspread[apv] += cal_->chi2_[apv][strip]*cal_->chi2_[apv][strip]/nStrips;
  }
                
  // fill the mean, max, min, spread, ... histograms.
  for(int i=0;i<2;++i) {
    cal_->mean_amplitude_[i] = Amean[i];
    cal_->mean_tail_[i] = Tmean[i];
    cal_->mean_riseTime_[i] = Rmean[i];
    cal_->mean_timeConstant_[i] = Cmean[i];
    cal_->mean_smearing_[i] = Smean[i];
    cal_->mean_chi2_[i] = Kmean[i];
    cal_->min_amplitude_[i] = Amin[i];
    cal_->min_tail_[i] = Tmin[i];
    cal_->min_riseTime_[i] = Rmin[i];
    cal_->min_timeConstant_[i] = Cmin[i];
    cal_->min_smearing_[i] = Smin[i];
    cal_->min_chi2_[i] = Kmin[i];
    cal_->max_amplitude_[i] = Amax[i];
    cal_->max_tail_[i] = Tmax[i];
    cal_->max_riseTime_[i] = Rmax[i];
    cal_->max_timeConstant_[i] = Cmax[i];
    cal_->max_smearing_[i] = Smax[i];
    cal_->max_chi2_[i] = Kmax[i];
    cal_->spread_amplitude_[i] = sqrt(fabs(Aspread[i]-Amean[i]*Amean[i]));
    cal_->spread_tail_[i] = sqrt(fabs(Tspread[i]-Tmean[i]*Tmean[i]));
    cal_->spread_riseTime_[i] = sqrt(fabs(Rspread[i]-Rmean[i]*Rmean[i]));
    cal_->spread_timeConstant_[i] = sqrt(fabs(Cspread[i]-Cmean[i]*Cmean[i]));
    cal_->spread_smearing_[i] = sqrt(fabs(Sspread[i]-Smean[i]*Smean[i]));
    cal_->spread_chi2_[i] = sqrt(fabs(Kspread[i]-Kmean[i]*Kmean[i]));
  }
}

void CalibrationAlgorithm::correctDistribution ( TH1 *  histo ) const

private

Definition at line 251 of file CalibrationAlgorithm.cc.

References cal_, and CalibrationAnalysis::isScan_.

Referenced by analyse(), and CalibrationAlgorithm().

{
  // return the curve
  histo->Scale(-1);
  if ( cal_ ) { if( cal_->isScan_ ) histo->Scale(1/16.); }
}

void CalibrationAlgorithm::extract ( const std::vector< TH1 * > &  )

overrideprivatevirtual

Extracts and organises histograms.

Implements CommissioningAlgorithm.

Definition at line 46 of file CalibrationAlgorithm.cc.

References CommissioningAnalysis::addErrorCode(), CommissioningAlgorithm::anal(), cal_, sistrip::CALIBRATION, sistrip::CALIBRATION_DECO, sistrip::CALIBRATION_SCAN, sistrip::CALIBRATION_SCAN_DECO, CommissioningAlgorithm::extractFedKey(), CommissioningAnalysis::fedKey(), histo_, CalibrationAnalysis::isScan_, sistrip::mlCommissioning_, sistrip::numberOfHistos_, fftjetcommon_cfi::title, tmp, and sistrip::unexpectedTask_.

Referenced by CalibrationAlgorithm().

                                                                 {
  if ( !anal() ) {
    edm::LogWarning(mlCommissioning_)
      << "[CalibrationAlgorithm::" << __func__ << "]"
      << " NULL pointer to base Analysis object!";
    return; 
  }

  CommissioningAnalysis* tmp = const_cast<CommissioningAnalysis*>( anal() );
  cal_ = dynamic_cast<CalibrationAnalysis*>( tmp );
  if ( !cal_ ) {
    edm::LogWarning(mlCommissioning_)
      << "[CalibrationAlgorithm::" << __func__ << "]"
      << " NULL pointer to derived Analysis object!";
    return; 
  }

  // Check number of histograms
  if ( histos.size() != 32 && histos.size() !=2 ) {
    cal_->addErrorCode(sistrip::numberOfHistos_);
  }
  
  // Extract FED key from histo title
  if ( !histos.empty() ) { cal_->fedKey( extractFedKey( histos.front() ) ); }
  
  // Extract histograms
  std::vector<TH1*>::const_iterator ihis = histos.begin();
  unsigned int cnt = 0;
  for ( ; ihis != histos.end(); ihis++,cnt++ ) {
    
    // Check for NULL pointer
    if ( !(*ihis) ) { continue; }
    
    // Check name
    SiStripHistoTitle title( (*ihis)->GetName() );
    if ( title.runType() != sistrip::CALIBRATION && 
     title.runType() != sistrip::CALIBRATION_DECO &&
         title.runType() != sistrip::CALIBRATION_SCAN && 
     title.runType() != sistrip::CALIBRATION_SCAN_DECO ) {
      cal_->addErrorCode(sistrip::unexpectedTask_);
      continue;
    }
    
    cal_->isScan_ = ( title.runType() == sistrip::CALIBRATION_SCAN || 
              title.runType() == sistrip::CALIBRATION_SCAN_DECO ); 
    
    // Extract calibration histo
    histo_[cnt].first = *ihis;
    histo_[cnt].second = (*ihis)->GetTitle();
  }
  
}

TF1 * CalibrationAlgorithm::fitPulse ( TH1 *  histo, float  rangeLow = 0, float  rangeHigh = -1  )

private

Definition at line 260 of file CalibrationAlgorithm.cc.

References cal_, CalibrationAnalysis::deconv_, deconv_fitter_, relativeConstraints::error, mps_fire::i, N, peak_fitter_, and mathSSE::sqrt().

Referenced by analyse(), CalibrationAlgorithm(), and maximum().

{
  if(!cal_) return nullptr;
  if(!histo) return nullptr;
  float noise = 4.;
  float N = round(histo->GetMaximum()/125.);
  float error = sqrt(2*N)*noise;
  //float error = sqrt(2)*noise*N; // ?????????
  for(int i=1;i<=histo->GetNbinsX();++i) {
    histo->SetBinError(i,error);
  }
  if (cal_->deconv_) {
    if(rangeLow>rangeHigh)
      histo->Fit(deconv_fitter_,"Q");
    else 
      histo->Fit(deconv_fitter_,"0Q","",rangeLow,rangeHigh);
    return deconv_fitter_; 
  } else {
    if(rangeLow>rangeHigh)
      histo->Fit(peak_fitter_,"Q");
    else 
      histo->Fit(peak_fitter_,"0Q","",rangeLow,rangeHigh);
    return peak_fitter_; 
  } 
}

const Histo& CalibrationAlgorithm::histo ( int  i ) const

inline

Definition at line 26 of file CalibrationAlgorithm.h.

References histo_, and mps_fire::i.

{ return histo_[i]; }

float CalibrationAlgorithm::maximum ( TH1 *  h )

private

Definition at line 290 of file CalibrationAlgorithm.cc.

References stringResolutionProvider_cfi::bin, and fitPulse().

Referenced by analyse(), and CalibrationAlgorithm().

                                            {
    int bin = h->GetMaximumBin();
    // fit around the maximum with the detector response and take the max from the fit
    TF1* fit = fitPulse(h,h->GetBinCenter(bin)-25,h->GetBinCenter(bin)+25);
    return fit->GetMaximumX();
}

float CalibrationAlgorithm::turnOn ( TH1 *  h )

private

Definition at line 299 of file CalibrationAlgorithm.cc.

References CustomPhysics_cfi::amplitude, runEdmFileComparison::base, stringResolutionProvider_cfi::bin, end, and ntuplemaker::time.

Referenced by analyse(), and CalibrationAlgorithm().

                                           {
    // localize the rising edge
    int bin=1;
    float amplitude = h->GetMaximum();
    for(; bin<= h->GetNbinsX() && h->GetBinContent(bin)<0.4*amplitude; ++bin) {}
    float end = h->GetBinLowEdge(bin);
    // fit the rising edge with a sigmoid
    TF1* sigmoid = new TF1("sigmoid","[0]/(1+exp(-[1]*(x-[2])))+[3]",0,end);
    sigmoid->SetParLimits(0,amplitude/10.,amplitude);
    sigmoid->SetParLimits(1,0.05,0.5);
    sigmoid->SetParLimits(2,end-10,end+10);
    sigmoid->SetParLimits(3,-amplitude/10.,amplitude/10.);
    sigmoid->SetParameters(amplitude/2.,0.1,end,0.);
    h->Fit(sigmoid,"0QR");
    // return the point where the fit = 3% signal.
    float time = 0.;
    float base = sigmoid->GetMinimum(0,end);
    for(;time<end && (sigmoid->Eval(time)-base)<0.03*(amplitude-base);time += 0.1) {}
    delete sigmoid;
    return time-0.05;
}

Member Data Documentation

CalibrationAnalysis* CalibrationAlgorithm::cal_

private

Definition at line 57 of file CalibrationAlgorithm.h.

Referenced by analyse(), correctDistribution(), extract(), and fitPulse().

TF1* CalibrationAlgorithm::deconv_fitter_

private

Fitter in deconvolution mode

Definition at line 52 of file CalibrationAlgorithm.h.

Referenced by CalibrationAlgorithm(), and fitPulse().

Histo CalibrationAlgorithm::histo_[32]

private

pulse shape

Definition at line 49 of file CalibrationAlgorithm.h.

Referenced by analyse(), extract(), and histo().

TF1* CalibrationAlgorithm::peak_fitter_

private

Fitter in peak mode

Definition at line 55 of file CalibrationAlgorithm.h.

Referenced by CalibrationAlgorithm(), and fitPulse().