+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * *
+ * Author: The ALICE Off-line Project. *
+ * Contributors are mentioned in the code where appropriate. *
+ * *
+ * Permission to use, copy, modify and distribute this software and its *
+ * documentation strictly for non-commercial purposes is hereby granted *
+ * without fee, provided that the above copyright notice appears in all *
+ * copies and that both the copyright notice and this permission notice *
+ * appear in the supporting documentation. The authors make no claims *
+ * about the suitability of this software for any purpose. It is *
+ * provided "as is" without express or implied warranty. *
+ **************************************************************************/
+
+/* $Id$ */
+
+
//-----------------------------------------------------------------------------
-// $Header$
+//
//
//
// Origin: Marian Ivanov, Uni. of Bratislava, ivanov@fmph.uniba.sk
// Declaration of class AliTPCRF1D
//
//-----------------------------------------------------------------------------
-#include "TMath.h"
+
+//
+
+#include <Riostream.h>
+#include <TCanvas.h>
+#include <TClass.h>
+#include <TF2.h>
+#include <TH1.h>
+#include <TMath.h>
+#include <TPad.h>
+#include <TString.h>
+#include <TStyle.h>
+
#include "AliTPCRF1D.h"
-#include "TF2.h"
-#include <iostream.h>
-#include "TCanvas.h"
-#include "TPad.h"
-#include "TStyle.h"
-#include "TH1.h"
-extern TStyle * gStyle;
+extern TStyle * gStyle;
+
+Int_t AliTPCRF1D::fgNRF=100; //default number of interpolation points
+Float_t AliTPCRF1D::fgRFDSTEP=0.01; //default step in cm
static Double_t funGauss(Double_t *x, Double_t * par)
{
+ //Gauss function -needde by the generic function object
return TMath::Exp(-(x[0]*x[0])/(2*par[0]*par[0]));
}
static Double_t funCosh(Double_t *x, Double_t * par)
{
+ //Cosh function -needde by the generic function object
return 1/TMath::CosH(3.14159*x[0]/(2*par[0]));
}
static Double_t funGati(Double_t *x, Double_t * par)
{
- //par[1] = is equal to k3
- //par[0] is equal to pad wire distance
- Float_t K3=par[1];
- Float_t K3R=TMath::Sqrt(K3);
- Float_t K2=(TMath::Pi()/2)*(1-K3R/2.);
- Float_t K1=K2*K3R/(4*TMath::ATan(K3R));
+ //Gati function -needde by the generic function object
+ Float_t k3=par[1];
+ Float_t k3R=TMath::Sqrt(k3);
+ Float_t k2=(TMath::Pi()/2)*(1-k3R/2.);
+ Float_t k1=k2*k3R/(4*TMath::ATan(k3R));
Float_t l=x[0]/par[0];
- Float_t tan2=TMath::TanH(K2*l);
+ Float_t tan2=TMath::TanH(k2*l);
tan2*=tan2;
- Float_t res = K1*(1-tan2)/(1+K3*tan2);
+ Float_t res = k1*(1-tan2)/(1+k3*tan2);
return res;
}
-
-
-
-///////////////////////////////////////////////////////////////////////////
-///////////////////////////////////////////////////////////////////////////
///////////////////////////////////////////////////////////////////////////
///////////////////////////////////////////////////////////////////////////
-AliTPCRF1D * gRF1D;
ClassImp(AliTPCRF1D)
AliTPCRF1D::AliTPCRF1D(Bool_t direct,Int_t np,Float_t step)
+ :TObject(),
+ fNRF(0),
+ fDSTEPM1(0.),
+ fcharge(0),
+ forigsigma(0.),
+ fpadWidth(3.5),
+ fkNorm(0.5),
+ fInteg(0.),
+ fGRF(0),
+ fSigma(0.),
+ fOffset(0.),
+ fDirect(kFALSE),
+ fPadDistance(0.)
{
+ //default constructor for response function object
fDirect=direct;
- fNRF = np;
+ if (np!=0)fNRF = np;
+ else (fNRF=fgNRF);
fcharge = new Float_t[fNRF];
- fDSTEPM1=1./step;
- fSigma = 0;
- gRF1D = this;
- fGRF = 0;
- fkNorm = 0.5;
- fpadWidth = 3.5;
- forigsigma=0.;
- fOffset = 0.;
+ if (step>0) fDSTEPM1=1./step;
+ else fDSTEPM1 = 1./fgRFDSTEP;
+ for(Int_t i=0;i<5;i++) {
+ funParam[i]=0.;
+ fType[i]=0;
+ }
+
}
+AliTPCRF1D::AliTPCRF1D(const AliTPCRF1D &prf)
+ :TObject(prf),
+ fNRF(prf.fNRF),
+ fDSTEPM1(prf.fDSTEPM1),
+ fcharge(0),
+ forigsigma(prf.forigsigma),
+ fpadWidth(prf.fpadWidth),
+ fkNorm(prf.fkNorm),
+ fInteg(prf.fInteg),
+ fGRF(new TF1(*(prf.fGRF))),
+ fSigma(prf.fSigma),
+ fOffset(prf.fOffset),
+ fDirect(prf.fDirect),
+ fPadDistance(prf.fPadDistance)
+{
+ //
+ //
+ for(Int_t i=0;i<5;i++) {
+ funParam[i]=0.;
+ fType[i]=0;
+ }
+ fcharge = new Float_t[fNRF];
+ memcpy(fcharge,prf.fcharge, fNRF*sizeof(Float_t));
+
+ //PH Change the name (add 0 to the end)
+ TString s(fGRF->GetName());
+ s+="0";
+ fGRF->SetName(s.Data());
+}
+
+AliTPCRF1D & AliTPCRF1D::operator = (const AliTPCRF1D &prf)
+{
+ if(this!=&prf) {
+ TObject::operator=(prf);
+ fNRF=prf.fNRF;
+ fDSTEPM1=prf.fDSTEPM1;
+ delete [] fcharge;
+ fcharge = new Float_t[fNRF];
+ memcpy(fcharge,prf.fcharge, fNRF*sizeof(Float_t));
+ forigsigma=prf.forigsigma;
+ fpadWidth=prf.fpadWidth;
+ fkNorm=prf.fkNorm;
+ fInteg=prf.fInteg;
+ delete fGRF;
+ fGRF=new TF1(*(prf.fGRF));
+ //PH Change the name (add 0 to the end)
+ TString s(fGRF->GetName());
+ s+="0";
+ fGRF->SetName(s.Data());
+ fSigma=prf.fSigma;
+ fOffset=prf.fOffset;
+ fDirect=prf.fDirect;
+ fPadDistance=prf.fPadDistance;
+ }
+ return *this;
+}
+
+
AliTPCRF1D::~AliTPCRF1D()
{
- if (fcharge!=0) delete fcharge;
- if (fGRF !=0 ) fGRF->Delete();
+ //
+ delete [] fcharge;
+ delete fGRF;
}
Float_t AliTPCRF1D::GetRF(Float_t xin)
{
- //x xin DSTEP unit
+ //function which return response
+ //for the charge in distance xin
//return linear aproximation of RF
- Float_t x = TMath::Abs((xin-fOffset)*fDSTEPM1)+fNRF/2;
+ Float_t x = (xin-fOffset)*fDSTEPM1+fNRF/2;
Int_t i1=Int_t(x);
if (x<0) i1-=1;
Float_t res=0;
- if (i1+1<fNRF)
+ if (i1+1<fNRF &&i1>0)
res = fcharge[i1]*(Float_t(i1+1)-x)+fcharge[i1+1]*(x-Float_t(i1));
return res;
}
Float_t AliTPCRF1D::GetGRF(Float_t xin)
{
+ //function which returnoriginal charge distribution
+ //this function is just normalised for fKnorm
if (fGRF != 0 )
return fkNorm*fGRF->Eval(xin)/fInteg;
else
void AliTPCRF1D::SetParam( TF1 * GRF,Float_t padwidth,
Float_t kNorm, Float_t sigma)
{
+ //adjust parameters of the original charge distribution
+ //and pad size parameters
fpadWidth = padwidth;
fGRF = GRF;
fkNorm = kNorm;
if (sigma==0) sigma= fpadWidth/TMath::Sqrt(12.);
forigsigma=sigma;
fDSTEPM1 = 10/TMath::Sqrt(sigma*sigma+fpadWidth*fpadWidth/12);
- sprintf(fType,"User");
+ //sprintf(fType,"User");
+ snprintf(fType,5,"User");
// Update();
}
void AliTPCRF1D::SetGauss(Float_t sigma, Float_t padWidth,
Float_t kNorm)
{
- // char s[120];
+ //
+ // set parameters for Gauss generic charge distribution
+ //
fpadWidth = padWidth;
fkNorm = kNorm;
if (fGRF !=0 ) fGRF->Delete();
- fGRF = new TF1("fun",funGauss,-5,5,2);
+ fGRF = new TF1("funGauss",funGauss,-5,5,1);
funParam[0]=sigma;
forigsigma=sigma;
fGRF->SetParameters(funParam);
fDSTEPM1 = 10./TMath::Sqrt(sigma*sigma+fpadWidth*fpadWidth/12);
- //by default I set the step as one tenth of sigma
- // Update();
- sprintf(fType,"Gauss");
+ //by default I set the step as one tenth of sigma
+ //sprintf(fType,"Gauss");
+ snprintf(fType,5,"Gauss");
}
void AliTPCRF1D::SetCosh(Float_t sigma, Float_t padWidth,
Float_t kNorm)
{
- // char s[120];
+ //
+ // set parameters for Cosh generic charge distribution
+ //
fpadWidth = padWidth;
fkNorm = kNorm;
if (fGRF !=0 ) fGRF->Delete();
- fGRF = new TF1("fun", funCosh, -5.,5.,2);
+ fGRF = new TF1("funCosh", funCosh, -5.,5.,2);
funParam[0]=sigma;
fGRF->SetParameters(funParam);
forigsigma=sigma;
fDSTEPM1 = 10./TMath::Sqrt(sigma*sigma+fpadWidth*fpadWidth/12);
//by default I set the step as one tenth of sigma
- // Update();
- sprintf(fType,"Cosh");
+ //sprintf(fType,"Cosh");
+ snprintf(fType,5,"Cosh");
}
void AliTPCRF1D::SetGati(Float_t K3, Float_t padDistance, Float_t padWidth,
Float_t kNorm)
{
- // char s[120];
+ //
+ // set parameters for Gati generic charge distribution
+ //
fpadWidth = padWidth;
fkNorm = kNorm;
if (fGRF !=0 ) fGRF->Delete();
- fGRF = new TF1("fun", funGati, -5.,5.,2);
+ fGRF = new TF1("funGati", funGati, -5.,5.,2);
funParam[0]=padDistance;
funParam[1]=K3;
fGRF->SetParameters(funParam);
forigsigma=padDistance;
fDSTEPM1 = 10./TMath::Sqrt(padDistance*padDistance+fpadWidth*fpadWidth/12);
//by default I set the step as one tenth of sigma
- // Update();
- sprintf(fType,"Gati");
+ //sprintf(fType,"Gati");
+ snprintf(fType,5,"Gati");
}
-void AliTPCRF1D::Draw(Float_t x1,Float_t x2,Int_t N)
+
+
+void AliTPCRF1D::DrawRF(Float_t x1,Float_t x2,Int_t N)
{
+ //
+ //Draw prf in selected region <x1,x2> with nuber of diviision = n
+ //
char s[100];
TCanvas * c1 = new TCanvas("canRF","Pad response function",700,900);
c1->cd();
TPad * pad2 = new TPad("pad2RF","",0.05,0.05,0.95,0.45,21);
pad2->Draw();
- sprintf(s,"RF response function for %1.2f cm pad width",
- fpadWidth);
+ //sprintf(s,"RF response function for %1.2f cm pad width",
+ // fpadWidth);
+ snprintf(s,60,"RF response function for %1.2f cm pad width",fpadWidth);
pad1->cd();
TH1F * hRFo = new TH1F("hRFo","Original charge distribution",N+1,x1,x2);
pad2->cd();
void AliTPCRF1D::Update()
{
- //initialize to 0
+ //
+ //update fields with interpolated values for
+ //PRF calculation
+
+ //at the begining initialize to 0
for (Int_t i =0; i<fNRF;i++) fcharge[i] = 0;
if ( fGRF == 0 ) return;
fInteg = fGRF->Integral(-5*forigsigma,5*forigsigma,funParam,0.00001);
void AliTPCRF1D::Streamer(TBuffer &R__b)
{
- // Stream an object of class AliTPC.
-
+ // Stream an object of class AliTPCRF1D.
if (R__b.IsReading()) {
- Version_t R__v = R__b.ReadVersion(); if (R__v) { }
- TObject::Streamer(R__b);
- //read pad parameters
- R__b >> fpadWidth;
- //read charge parameters
- R__b >> fType[0];
- R__b >> fType[1];
- R__b >> fType[2];
- R__b >> fType[3];
- R__b >> fType[4];
- R__b >> forigsigma;
- R__b >> fkNorm;
- R__b >> fK3X;
- R__b >> fPadDistance;
- R__b >> fInteg;
- R__b >> fOffset;
+ AliTPCRF1D::Class()->ReadBuffer(R__b, this);
//read functions
- if (fGRF!=0) {
- delete fGRF;
- fGRF=0;
- }
- if (strncmp(fType,"User",3)==0){
- fGRF= new TF1;
- R__b>>fGRF;
- }
- if (strncmp(fType,"Gauss",3)==0)
- fGRF = new TF1("fun",funGauss,-5.,5.,4);
- if (strncmp(fType,"Cosh",3)==0)
- fGRF = new TF1("fun",funCosh,-5.,5.,4);
- if (strncmp(fType,"Gati",3)==0)
- fGRF = new TF1("fun",funGati,-5.,5.,4);
- R__b >>fDSTEPM1;
- R__b >>fNRF;
- R__b.ReadFastArray(fcharge,fNRF);
- R__b.ReadFastArray(funParam,5);
- if (fGRF!=0) fGRF->SetParameters(funParam);
+ if (strncmp(fType,"Gauss",3)==0) {delete fGRF; fGRF = new TF1("funGauss",funGauss,-5.,5.,4);}
+ if (strncmp(fType,"Cosh",3)==0) {delete fGRF; fGRF = new TF1("funCosh",funCosh,-5.,5.,4);}
+ if (strncmp(fType,"Gati",3)==0) {delete fGRF; fGRF = new TF1("funGati",funGati,-5.,5.,4);}
+ if (fGRF) fGRF->SetParameters(funParam);
} else {
- R__b.WriteVersion(AliTPCRF1D::IsA());
- TObject::Streamer(R__b);
- //write pad width
- R__b << fpadWidth;
- //write charge parameters
- R__b << fType[0];
- R__b << fType[1];
- R__b << fType[2];
- R__b << fType[3];
- R__b << fType[4];
- R__b << forigsigma;
- R__b << fkNorm;
- R__b << fK3X;
- R__b << fPadDistance;
- R__b << fInteg;
- R__b << fOffset;
- //write interpolation parameters
- if (strncmp(fType,"User",3)==0) R__b <<fGRF;
- R__b <<fDSTEPM1;
- R__b <<fNRF;
- R__b.WriteFastArray(fcharge,fNRF);
- R__b.WriteFastArray(funParam,5);
-
-
-
+ AliTPCRF1D::Class()->WriteBuffer(R__b, this);
}
}
+
+
+Double_t AliTPCRF1D::Gamma4(Double_t x, Double_t p0, Double_t p1){
+ //
+ // Gamma 4 Time response function of ALTRO
+ //
+ if (x<0) return 0;
+ Double_t g1 = TMath::Exp(-4.*x/p1);
+ Double_t g2 = TMath::Power(x/p1,4);
+ return p0*g1*g2;
+}