PedsFullNoiseAlgorithm Class Reference

Histogram-based analysis for pedestal run. More...

#include <PedsFullNoiseAlgorithm.h>

Inheritance diagram for PedsFullNoiseAlgorithm:

Public Member Functions
const Histo &	hNoise () const
const Histo &	hNoise1D () const
const Histo &	hPeds () const
	PedsFullNoiseAlgorithm (const edm::ParameterSet &pset, PedsFullNoiseAnalysis *const )
virtual	~PedsFullNoiseAlgorithm ()
Public Member Functions inherited from CommissioningAlgorithm
void	analysis (const std::vector< TH1 * > &)
	CommissioningAlgorithm (CommissioningAnalysis *const )
	CommissioningAlgorithm ()
virtual	~CommissioningAlgorithm ()

Private Member Functions
void	analyse ()
void	extract (const std::vector< TH1 * > &)
	PedsFullNoiseAlgorithm ()

Private Attributes
float	deadStripMax_
Histo	hNoise1D_
Histo	hNoise_
Histo	hPeds_
float	ksProbCut_
std::string	noiseDef_
float	noisyStripMin_

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

Histogram-based analysis for pedestal run.

Author

M. Wingham, R.Bainbridge

Definition at line 16 of file PedsFullNoiseAlgorithm.h.

Constructor & Destructor Documentation

PedsFullNoiseAlgorithm::PedsFullNoiseAlgorithm	(	const edm::ParameterSet &	pset,
		PedsFullNoiseAnalysis * const	anal
	)

Definition at line 18 of file PedsFullNoiseAlgorithm.cc.

  : CommissioningAlgorithm(anal),
    hPeds_(0,""),
    hNoise_(0,""),
    hNoise1D_(0,""),
    deadStripMax_(pset.getParameter<double>("DeadStripMax")),
    noisyStripMin_(pset.getParameter<double>("NoisyStripMin")),
    noiseDef_(pset.getParameter<std::string>("NoiseDefinition")),
    ksProbCut_(pset.getParameter<double>("KsProbCut"))
{
  //LogDebug(mlCommissioning_)
  //  << "[PedsFullNoiseAlgorithm::" << __func__ << "]"
  //  << " Set maximum noise deviation for dead strip determination to: " << deadStripMax_;
    // LogDebug(mlCommissioning_)
  //  << "[PedsFullNoiseAlgorithm::" << __func__ << "]"
  //  << " Set minimal noise deviation for noise strip determination to: " << noisyStripMin_;
}

virtual PedsFullNoiseAlgorithm::~PedsFullNoiseAlgorithm ( )

inlinevirtual

Definition at line 22 of file PedsFullNoiseAlgorithm.h.

References hNoise(), hNoise1D(), and hPeds().

{;}

PedsFullNoiseAlgorithm::PedsFullNoiseAlgorithm ( )

inlineprivate

Definition at line 32 of file PedsFullNoiseAlgorithm.h.

References analyse(), and extract().

{;}

Member Function Documentation

void PedsFullNoiseAlgorithm::analyse ( )

privatevirtual

Performs histogram anaysis.

Implements CommissioningAlgorithm.

Definition at line 97 of file PedsFullNoiseAlgorithm.cc.

References a, CommissioningAnalysis::addErrorCode(), CommissioningAlgorithm::anal(), b, createTree::dd, PedsFullNoiseAnalysis::dead_, deadStripMax_, f, CommissioningAnalysis::fecKey(), hNoise_, hPeds_, mps_fire::i, reco::if(), sistrip::invalid_, PedsFullNoiseAnalysis::ksProb_, ksProbCut_, LogTrace, sistrip::maximum_, sistrip::mlCommissioning_, sistrip::mlDqmClient_, PedsFullNoiseAnalysis::noise_, PedsFullNoiseAnalysis::noiseBin84_, noiseDef_, PedsFullNoiseAnalysis::noiseGaus_, PedsFullNoiseAnalysis::noiseMax_, PedsFullNoiseAnalysis::noiseMean_, PedsFullNoiseAnalysis::noiseMin_, PedsFullNoiseAnalysis::noiseRMS_, PedsFullNoiseAnalysis::noiseSignif_, PedsFullNoiseAnalysis::noiseSpread_, PedsFullNoiseAnalysis::noisy_, noisyStripMin_, sistrip::nullPtr_, sistrip::numberOfBins_, PedsFullNoiseAnalysis::peds_, PedsFullNoiseAnalysis::pedsMax_, PedsFullNoiseAnalysis::pedsMean_, PedsFullNoiseAnalysis::pedsMin_, PedsFullNoiseAnalysis::pedsSpread_, PedsFullNoiseAnalysis::raw_, PedsFullNoiseAnalysis::rawMax_, PedsFullNoiseAnalysis::rawMean_, PedsFullNoiseAnalysis::rawMin_, PedsFullNoiseAnalysis::rawSpread_, mathSSE::sqrt(), tmp, w, and x.

Referenced by PedsFullNoiseAlgorithm().

                                     {

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

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

  if ( !hPeds_.first ) {
    ana->addErrorCode(sistrip::nullPtr_);
    return;
  }

  if ( !hNoise_.first ) {
    ana->addErrorCode(sistrip::nullPtr_);
    return;
  }
  
  TProfile * peds_histo = dynamic_cast<TProfile *>(hPeds_.first);
  TH2S * noise_histo        = dynamic_cast<TH2S *>(hNoise_.first);
  if ( !peds_histo ) {
    ana->addErrorCode(sistrip::nullPtr_);
    return;
  }

  if ( !noise_histo ) {
    ana->addErrorCode(sistrip::nullPtr_);
    return;
  }

  if ( peds_histo->GetNbinsX() != 256 ) {
    ana->addErrorCode(sistrip::numberOfBins_);
    return;
  }

  if ( noise_histo->GetNbinsY() != 256 ) { // X range is configurable
    ana->addErrorCode(sistrip::numberOfBins_);
    return;
  }
  
    // Iterate through APVs 
    for ( uint16_t iapv = 0; iapv < 2; iapv++ ) {
    
    // Used to calc mean and rms for peds and noise
    float p_sum = 0., p_sum2 = 0., p_max = -1.*sistrip::invalid_, p_min = sistrip::invalid_;
    float n_sum = 0., n_sum2 = 0., n_max = -1.*sistrip::invalid_, n_min = sistrip::invalid_;
    float r_sum = 0., r_sum2 = 0., r_max = -1.*sistrip::invalid_, r_min = sistrip::invalid_;           
    // Iterate through strips of APV

    //Need our own copy for thread safety
    TF1 mygaus("mygaus","gaus");
    for ( uint16_t istr = 0; istr < 128; istr++ ) {
        
      ana->ksProb_[iapv].push_back(0);
      ana->noiseGaus_[iapv].push_back(0);
      ana->noiseBin84_[iapv].push_back(0);
      ana->noiseRMS_[iapv].push_back(0);
      ana->noiseSignif_[iapv].push_back(0);
    
      // pedestals and raw noise
      if ( peds_histo ) {
            if ( peds_histo->GetBinEntries(iapv*128 + istr + 1) ) {
                ana->peds_[iapv][istr] = peds_histo->GetBinContent(iapv*128 + istr + 1);
                p_sum += ana->peds_[iapv][istr];
                p_sum2 += (ana->peds_[iapv][istr] * ana->peds_[iapv][istr]);
                if ( ana->peds_[iapv][istr] > p_max ) { p_max = ana->peds_[iapv][istr];}
                if ( ana->peds_[iapv][istr] < p_min ) { p_min = ana->peds_[iapv][istr];}
                ana->raw_[iapv][istr] = peds_histo->GetBinError(iapv*128 + istr + 1);
                r_sum += ana->raw_[iapv][istr];
                r_sum2 += (ana->raw_[iapv][istr] * ana->raw_[iapv][istr]);
                if ( ana->raw_[iapv][istr] > r_max ) { r_max = ana->raw_[iapv][istr]; }
                if ( ana->raw_[iapv][istr] < r_min ) { r_min = ana->raw_[iapv][istr]; }
            }
      }
      // Noise from Full Distribution
      if ( noise_histo ) {
        // Fit the ADC Distribution from TH2S by projecting it out and fitting.
        
        TH1S * noisehist = new TH1S("noisehist","",noise_histo->GetNbinsX(),
                                -noise_histo->GetNbinsX()/2,noise_histo->GetNbinsX()/2);
                    
        
        for(int i=0;i<=noise_histo->GetNbinsX()+1;i++){
         noisehist->SetBinContent(i,noise_histo->GetBinContent(i,iapv*128 + istr + 1));
        }           
        // If the histogram is valid.
        if(noisehist->Integral() > 0){
            ana->noiseRMS_[iapv][istr] = noisehist->GetRMS();
            noisehist->Fit(&mygaus,"Q");                       
            ana->noiseGaus_[iapv][istr]     = mygaus.GetParameter(2);
          
          // new Bin84 method
            std::vector<float> integralFrac;
                integralFrac.push_back(1.*noisehist->GetBinContent(0)/noisehist->Integral(0,noisehist->GetNbinsX()));
            // Calculate the integral of distro as a function of bins.
            for(int i = 1; i < noisehist->GetNbinsX();i++){
                integralFrac.push_back(float(noisehist->GetBinContent(i))/
                noisehist->Integral(0,noisehist->GetNbinsX())+integralFrac[i-1]);
                //Take the two bins next to 84% and solve for x in 0.84 = mx+b
                if (integralFrac[i] >= 0.84135 && integralFrac[i-1] < 0.84135) {
                // my quadratic noise calculation
                    double w = noisehist->GetBinWidth(i);
                    double a = noisehist->GetBinContent(i-1);
                    double b = noisehist->GetBinContent(i);
                    double f = w*(0.84135 -integralFrac[i-1])/(integralFrac[i]-integralFrac[i-1]);
                    double x = 0;
                    if (a==b) {
                        x = f;
                    } else {
                        double aa = (b-a)/(2*w);
                        double bb = (b+a)/2;
                        double cc = -b*f;
                        double dd = bb*bb-4*aa*cc; //if (dd<0) dd=0;
                        x = (-bb+sqrt(dd))/(2*aa);
                    }
                    ana->noiseBin84_[iapv][istr] = noisehist->GetBinLowEdge(i) + x;
                }
          } 
            // Compare shape of ADC to a histogram made of Gaus Fit for KSProb, Chi2Prob, Etc...
            TH1D * FitHisto = new TH1D("FitHisto","FitHisto",noisehist->GetNbinsX(),
                           -noisehist->GetNbinsX()/2,noisehist->GetNbinsX()/2);
            FitHisto->Add(&mygaus);
            FitHisto->Sumw2();
            noisehist->Sumw2();
          
            if(FitHisto->Integral() > 0){
                // This is stored as a float but will be plotted as an int at the summary histos.
            // This forces the histo to draw 10000 bins instead of 1.
                ana->ksProb_[iapv][istr] = noisehist->KolmogorovTest(FitHisto)*10000;
            }               
            delete FitHisto;
        }
        delete noisehist;    
      }    
      // Assigning the actual noise values used for Upload!!!!!!!!!!!!!!!!!!!!
      if (noiseDef_ == "Bin84") {
        if (ana->noiseBin84_[iapv][istr] > 0) {
          ana->noise_[iapv][istr] = ana->noiseBin84_[iapv][istr];
        } else {
          ana->noise_[iapv][istr] = ana->noiseRMS_[iapv][istr];
        }
      } else if (noiseDef_ == "RMS") {
        ana->noise_[iapv][istr] = ana->noiseRMS_[iapv][istr];
      } else edm::LogWarning(mlCommissioning_)<< "[PedsFullNoiseAlgorithm::"
                    << __func__ << "]"<< " Unknown noise definition!!!";
      
      // Setting Sum of RMS and RMS^2 for Dead/Noisy Strip calculations
      n_sum += ana->noise_[iapv][istr];
      n_sum2 += (ana->noise_[iapv][istr] * ana->noise_[iapv][istr]);
      if ( ana->noise_[iapv][istr] > n_max ) { n_max = ana->noise_[iapv][istr]; }
      if ( ana->noise_[iapv][istr] < n_min ) { n_min = ana->noise_[iapv][istr]; }
      
    } // strip loop
  
    // Calc mean and rms for peds
    if ( !ana->peds_[iapv].empty() ) { 
        p_sum /= static_cast<float>( ana->peds_[iapv].size() );
      p_sum2 /= static_cast<float>( ana->peds_[iapv].size() );
      ana->pedsMean_[iapv] = p_sum;
      ana->pedsSpread_[iapv] = sqrt( fabs(p_sum2 - p_sum*p_sum) );
    }
  
    // Calc mean and rms for noise using noiseRMS.
    if ( !ana->noise_[iapv].empty() ) { 
        n_sum /= static_cast<float>( ana->noise_[iapv].size() );
      n_sum2 /= static_cast<float>( ana->noise_[iapv].size() );
      ana->noiseMean_[iapv] = n_sum;
      ana->noiseSpread_[iapv] = sqrt( fabs(n_sum2 - n_sum*n_sum) );
    }

    // Calc mean and rms for raw noise
    if ( !ana->raw_[iapv].empty() ) { 
        r_sum /= static_cast<float>( ana->raw_[iapv].size() );
      r_sum2 /= static_cast<float>( ana->raw_[iapv].size() );
      ana->rawMean_[iapv] = r_sum;
      ana->rawSpread_[iapv] = sqrt( fabs(r_sum2 - r_sum*r_sum) );
    }
      
    // Set max and min values for peds, noise and raw noise
    if ( p_max > -1.*sistrip::maximum_ ) { ana->pedsMax_[iapv] = p_max; }
    if ( p_min < 1.*sistrip::maximum_ )  { ana->pedsMin_[iapv] = p_min; }
    if ( n_max > -1.*sistrip::maximum_ ) { ana->noiseMax_[iapv] = n_max; }
    if ( n_min < 1.*sistrip::maximum_ )  { ana->noiseMin_[iapv] = n_min; }
    if ( r_max > -1.*sistrip::maximum_ ) { ana->rawMax_[iapv] = r_max; }
    if ( r_min < 1.*sistrip::maximum_ )  { ana->rawMin_[iapv] = r_min; }
  
    // Set dead and noisy strips
    for ( uint16_t istr = 0; istr < 128; istr++ ) { // strip loop 
      // Set the significance of the noise of each strip also compared to apv avg.
      ana->noiseSignif_[iapv][istr] = (ana->noise_[iapv][istr]-ana->noiseMean_[iapv])/ana->noiseSpread_[iapv];
      if ( ana->noiseMin_[iapv] > sistrip::maximum_ || ana->noiseMax_[iapv] > sistrip::maximum_ ) { 
        continue; 
      }
      // Strip Masking for Dead Strips
      if(ana->noiseSignif_[iapv][istr] < -deadStripMax_){
        ana->dead_[iapv].push_back(istr);
        SiStripFecKey fec_key(ana->fecKey());
        LogTrace(mlDqmClient_)<<"DeadSignif "<<ana->noiseSignif_[iapv][istr]
        <<" "<<fec_key.fecCrate()
        <<" "<<fec_key.fecSlot()
        <<" "<<fec_key.fecRing()
        <<" "<<fec_key.ccuAddr()
        <<" "<<fec_key.ccuChan()
        <<" "<<fec_key.lldChan()
        <<" "<<iapv*128+istr<<std::endl;
      } // Strip Masking for Dead Strips
      // Laurent's Method for Flagging bad strips in TIB
      else if((ana->noiseMax_[iapv]/ana->noiseMean_[iapv] > 3 || ana->noiseSpread_[iapv] > 3)
      && ana->noiseSignif_[iapv][istr] > 1){
        ana->noisy_[iapv].push_back(istr);
        SiStripFecKey fec_key(ana->fecKey());
        LogTrace(mlDqmClient_)<<"NoisyLM "<<ana->noiseMax_[iapv]/ana->noiseMean_[iapv]     
        <<" "<<fec_key.fecCrate()
        <<" "<<fec_key.fecSlot()
        <<" "<<fec_key.fecRing()
        <<" "<<fec_key.ccuAddr()
        <<" "<<fec_key.ccuChan()
        <<" "<<fec_key.lldChan()
        <<" "<<iapv*128+istr<<std::endl;
      } // if NoisyLM 
      //Strip Masking for Non Gassian Strips which are also noisy.
      else if(ana->ksProb_[iapv][istr] < ksProbCut_){
        ana->noisy_[iapv].push_back(istr);
        SiStripFecKey fec_key(ana->fecKey());  
        LogTrace(mlDqmClient_)<<"NoisyKS "<<ana->ksProb_[iapv][istr] 
        <<" "<<fec_key.fecCrate()
        <<" "<<fec_key.fecSlot()
        <<" "<<fec_key.fecRing()
        <<" "<<fec_key.ccuAddr()
        <<" "<<fec_key.ccuChan()
        <<" "<<fec_key.lldChan()
        <<" "<<iapv*128+istr<<std::endl;
      } //Strip Masking for Non Gassian Strips which are also noisy.
      else if(ana->noiseSignif_[iapv][istr] > noisyStripMin_){
        ana->noisy_[iapv].push_back(istr);
        SiStripFecKey fec_key(ana->fecKey());
        LogTrace(mlDqmClient_)<<"NoisySignif "<<ana->noiseSignif_[iapv][istr]    
        <<" "<<fec_key.fecCrate()
        <<" "<<fec_key.fecSlot()
        <<" "<<fec_key.fecRing()
        <<" "<<fec_key.ccuAddr()
        <<" "<<fec_key.ccuChan()
        <<" "<<fec_key.lldChan()
        <<" "<<iapv*128+istr<<std::endl;
      } // if Signif 
    }// strip loop to set dead or noisy strips   
    } // apv loop
  //std::cout << std::endl;
}

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

privatevirtual

Extracts and organises histograms.

Implements CommissioningAlgorithm.

Definition at line 39 of file PedsFullNoiseAlgorithm.cc.

References CommissioningAnalysis::addErrorCode(), CommissioningAlgorithm::anal(), sistrip::extrainfo::commonMode_, CommissioningAlgorithm::extractFedKey(), CommissioningAnalysis::fedKey(), hNoise1D_, hNoise_, hPeds_, sistrip::mlCommissioning_, sistrip::extrainfo::noise2D_, sistrip::extrainfo::noiseProfile_, sistrip::numberOfHistos_, sistrip::extrainfo::pedestals_, sistrip::extrainfo::roughPedestals_, fftjetcommon_cfi::title, and sistrip::unexpectedExtraInfo_.

Referenced by PedsFullNoiseAlgorithm().

                                                                    { 

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

  // Check number of histograms
  if ( histos.size() != 3 ) { 
    anal()->addErrorCode(sistrip::numberOfHistos_);
  }
  
  // Extract FED key from histo title
  if ( !histos.empty() ) { 
    anal()->fedKey( extractFedKey( histos.front() ) );
  }
  
  // Extract 1D histograms
  std::vector<TH1*>::const_iterator ihis = histos.begin();
  for ( ; ihis != histos.end(); ihis++ ) {
    
    // Check for NULL pointer
    if ( !(*ihis) ) { continue; }
    
// TO BE UPDATED!!!
    // Check run type
    SiStripHistoTitle title( (*ihis)->GetName() );
/*    if ( title.runType() != sistrip::PEDS_FULL_NOISE ) {
      anal()->addErrorCode(sistrip::unexpectedTask_);
      continue;
    }
*/
    // Extract peds histos
    if ( title.extraInfo().find(sistrip::extrainfo::roughPedestals_) != std::string::npos ) {
      //@@ something here for rough peds?
    } else if ( title.extraInfo().find(sistrip::extrainfo::pedestals_) != std::string::npos ) {
      hPeds_.first = *ihis;
      hPeds_.second = (*ihis)->GetName();
    } else if ( title.extraInfo().find(sistrip::extrainfo::commonMode_) != std::string::npos ) {
      //@@ something here for CM plots?
    } else if ( title.extraInfo().find(sistrip::extrainfo::noiseProfile_) != std::string::npos ) {
      //@@ something here for noise profile plot?
      hNoise1D_.first = *ihis;
      hNoise1D_.second = (*ihis)->GetName();
    } else if ( title.extraInfo().find(sistrip::extrainfo::noise2D_) != std::string::npos ) {
      hNoise_.first = *ihis;
      hNoise_.second = (*ihis)->GetName();
    } else { 
      anal()->addErrorCode(sistrip::unexpectedExtraInfo_);
    }  
  }

}

const PedsFullNoiseAlgorithm::Histo & PedsFullNoiseAlgorithm::hNoise ( ) const

inline

Definition at line 58 of file PedsFullNoiseAlgorithm.h.

References hNoise_.

Referenced by ~PedsFullNoiseAlgorithm().

{ return hNoise_; }

const Histo& PedsFullNoiseAlgorithm::hNoise1D ( ) const

inline

Referenced by ~PedsFullNoiseAlgorithm().

const PedsFullNoiseAlgorithm::Histo & PedsFullNoiseAlgorithm::hPeds ( ) const

inline

Definition at line 56 of file PedsFullNoiseAlgorithm.h.

References hPeds_.

Referenced by ~PedsFullNoiseAlgorithm().

{ return hPeds_; }

Member Data Documentation

float PedsFullNoiseAlgorithm::deadStripMax_

private

Analysis parameters

Definition at line 50 of file PedsFullNoiseAlgorithm.h.

Referenced by analyse().

Histo PedsFullNoiseAlgorithm::hNoise1D_

private

Definition at line 47 of file PedsFullNoiseAlgorithm.h.

Referenced by extract().

Histo PedsFullNoiseAlgorithm::hNoise_

private

Residuals and noise

Definition at line 46 of file PedsFullNoiseAlgorithm.h.

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

Histo PedsFullNoiseAlgorithm::hPeds_

private

Pedestals and raw noise

Definition at line 43 of file PedsFullNoiseAlgorithm.h.

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

float PedsFullNoiseAlgorithm::ksProbCut_

private

Definition at line 53 of file PedsFullNoiseAlgorithm.h.

Referenced by analyse().

std::string PedsFullNoiseAlgorithm::noiseDef_

private

Definition at line 52 of file PedsFullNoiseAlgorithm.h.

Referenced by analyse().

float PedsFullNoiseAlgorithm::noisyStripMin_

private

Definition at line 51 of file PedsFullNoiseAlgorithm.h.

Referenced by analyse().