LArPhysWaveHECTool.cxx File Reference

#include "LArCalibUtils/LArPhysWaveHECTool.h"
#include "LArIdentifier/LArOnline_SuperCellID.h"
#include <gsl/gsl_integration.h>
#include <gsl/gsl_errno.h>
#include <TGraphSmooth.h>
#include <iostream>
#include <fstream>
#include <vector>
#include <cmath>

Include dependency graph for LArPhysWaveHECTool.cxx:

Go to the source code of this file.

Classes
struct	gsl_f_params
struct	usertim

Functions
Double_t	normalize (TF1 *func, Double_t *rampl=NULL, Double_t from=0., Double_t to=0., Double_t step=1.)
Double_t	normalize_prof (TProfile *histo, Double_t *rampl=NULL)
void	Smooth_new (Double_t *x, Double_t *y, Double_t *yref, Int_t ymax)
void	SmootKern (Double_t *x, Double_t *y, Int_t ymax, Double_t band=10.)

Variables
	ATLAS_NO_CHECK_FILE_THREAD_SAFETY
const int	NMAX =1000
int	LArWaveFlag =0
int	Tp4_gsl_err_counter =0
int	Tp5_gsl_err_counter =0

Function Documentation

normalize()

Double_t normalize	(	TF1 *	func,
		Double_t *	rampl = NULL,
		Double_t	from = 0.,
		Double_t	to = 0.,
		Double_t	step = 1.
	)

Definition at line 825 of file LArPhysWaveHECTool.cxx.

{
 static Double_t amax, tmax,t, ampl;
 static Double_t t1,t2,tt1,tt2;
 static Int_t i;
 
 if((from == to) || (from > to)) 
     func->GetRange(t1,t2);
 else {
    func->GetRange(tt1,tt2);
    if(from>tt1) t1 = from; else t1 = tt1;
    if(to<tt2) t2 = to; else t2 = tt2;
  }
 
 amax = 0; tmax = t1;
 for(t=t1; t<=t2; t+=step)
  {
   if(func->Eval(t) > amax) {
     amax = func->Eval(t); tmax = t; 
    }
   }
 Double_t x[20], y[20];
 i=0;
 for(t=tmax-5; t<tmax+5; ++t) {
    x[i] = t;
    y[i] = func->Eval(t);
    ++i;
   }
 TGraph *g = new TGraph(10,x,y);
 g->Fit("pol2","QIE0");
 tmax = - g->GetFunction("pol2")->GetParameter(1) / (2 *g->GetFunction("pol2")->GetParameter(2));
 amax = g->GetFunction("pol2")->Eval(tmax);
 
 ampl = func->GetParameter(1);
 if(amax != 0) func->SetParameter(1,ampl/amax);
   else tmax = 0;
 cout << "Normalize: tmax="<<tmax << " ampl=" << ampl<<" amax="<<amax << endl;
 if(rampl != NULL) *rampl = amax;
 g->Delete();
 return tmax;
}

normalize_prof()

Double_t normalize_prof	(	TProfile *	histo,
		Double_t *	rampl = NULL
	)

Definition at line 868 of file LArPhysWaveHECTool.cxx.

{
 Double_t amax, tmax, x, y, err;
 Int_t nbins,i, j, nent;
 amax = histo->GetBinContent(0);
 tmax = histo->GetBinCenter(0);
 nbins = histo->GetNbinsX();
 for(i=1; i<nbins; i++)
  {
   if(amax < histo->GetBinContent(i))
     {
      amax = histo->GetBinContent(i);
      tmax =  histo->GetBinCenter(i);
     }
  }
 
 histo->Fit("pol2","0QI","",tmax-10., tmax+10.);
 tmax = - histo->GetFunction("pol2")->GetParameter(1) / (2 * histo->GetFunction("pol2")->GetParameter(2));
 amax = histo->GetFunction("pol2")->Eval(tmax);
 TProfile *htmp = new TProfile(*histo);
 histo->Reset();
 for(i=0; i<=nbins; i++)
  {
   if(htmp->GetBinEntries(i) == 0) continue;
   x = htmp->GetBinCenter(i);
   y = htmp->GetBinContent(i);
   nent=int(htmp->GetBinEntries(i));
   err = htmp->GetBinError(i);
// Hack for stupid root v.3, which do not like the one entries in profile histo
   for(j=0; j<7; ++j) {
      if(nent != 0) {
           histo->Fill(x,y/amax+sqrt(2)*err/nent);
           histo->Fill(x,y/amax-sqrt(2)*err/nent);
      } else {
           histo->Fill(x,y/amax+sqrt(2)*err);
           histo->Fill(x,y/amax-sqrt(2)*err);
      }
    }
  }
 
 if(rampl != NULL) *rampl = amax;
 
 return tmax;
}

Smooth_new()

void Smooth_new	(	Double_t *	x,
		Double_t *	y,
		Double_t *	yref,
		Int_t	ymax
	)

Definition at line 925 of file LArPhysWaveHECTool.cxx.

{
// Smoothing the difference between the yref and y, 
// then y = y + smoothed differnce
// repeating few times
  static Double_t *dif = NULL;
  static Int_t i,j;
  
  if(dif) delete[] dif;
  dif = new Double_t[ymax];
  for(j=0; j<4; ++j) {
    for(i=0; i<ymax; ++i) dif[i] = y[i] - yref[i];
    SmootKern(x,dif,ymax, 10.);
    for(i=0; i<ymax; ++i) yref[i] += dif[i];
  }
}

SmootKern()

void SmootKern	(	Double_t *	x,
		Double_t *	y,
		Int_t	ymax,
		Double_t	band = 10.
	)

Definition at line 912 of file LArPhysWaveHECTool.cxx.

{
   Double_t *xx, *yy;
   TGraph *g=new TGraph(ymax,x,y);
   TGraphSmooth *gs = new TGraphSmooth();
   TGraph *g1 = gs->SmoothKern(g,"normal",band);
   xx = g1->GetX(); yy = g1->GetY();
   for(Int_t i=0; i<g->GetN(); ++i) {
      x[i] = xx[i]; y[i] = yy[i];
   }
   delete g; delete gs;
}

Variable Documentation

ATLAS_NO_CHECK_FILE_THREAD_SAFETY

ATLAS_NO_CHECK_FILE_THREAD_SAFETY

Definition at line 17 of file LArPhysWaveHECTool.cxx.

LArWaveFlag

int LArWaveFlag =0

Definition at line 57 of file LArPhysWaveHECTool.cxx.

NMAX

const int NMAX =1000

Definition at line 22 of file LArPhysWaveHECTool.cxx.

Tp4_gsl_err_counter

int Tp4_gsl_err_counter =0

Definition at line 57 of file LArPhysWaveHECTool.cxx.

Tp5_gsl_err_counter

int Tp5_gsl_err_counter =0

Definition at line 57 of file LArPhysWaveHECTool.cxx.