Offline/HLT Offline/HLT OCDB entries (AliTPCClusterParam)
*/
+/*
+
+How to retrive it from file (created using calibration task):
+
+gSystem->Load("libANALYSIS");
+gSystem->Load("libTPCcalib");
+TFile fcalib("CalibObjects.root");
+TObjArray * array = (TObjArray*)fcalib.Get("TPCCalib");
+AliTPCcalibTracks * calibTracks = ( AliTPCcalibTracks *)array->FindObject("calibTracks");
+
+
+//USAGE of debug stream example
+ gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
+ gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
+ AliXRDPROOFtoolkit tool;
+ TChain * chainres = tool.MakeChain("tracks.txt","ResolCl",0,10200);
+ chainres->Lookup();
+*/
+
+
//
///////////////////////////////////////////////////////////////////////////////
#include <TChain.h>
#include <TFile.h>
#include <TH3F.h>
+#include <TProfile.h>
+
//
//#include <TPDGCode.h>
#include <TStyle.h>
#include <TCollection.h>
#include <iostream>
#include <TLinearFitter.h>
+#include <TString.h>
//
// AliROOT includes
#include "TText.h"
#include "TPaveText.h"
#include "TSystem.h"
+#include "TStatToolkit.h"
+#include "TCut.h"
-// Thread-stuff
-//#include "TThread.h"
-//#include "TMutex.h"
-//#include "TLockFile.h"
ClassImp(AliTPCcalibTracks)
AliTPCcalibTracks::AliTPCcalibTracks():
- TNamed(),
+ AliTPCcalibBase(),
fClusterParam(0),
- fDebugStream(0),
fROC(0),
fArrayAmpRow(0),
fArrayAmp(0),
fHclusterPerPadrowRaw(0),
fClusterCutHisto(0),
fCalPadClusterPerPad(0),
- fCalPadClusterPerPadRaw(0),
- fDebugLevel(0),
- fFitterLinY1(0), //!
- fFitterLinZ1(0), //!
- fFitterLinY2(0), //!
- fFitterLinZ2(0), //!
- fFitterParY(0), //!
- fFitterParZ(0) //!
+ fCalPadClusterPerPadRaw(0)
{
//
// AliTPCcalibTracks default constructor
//
- if (fDebugLevel > 0) cout << "AliTPCcalibTracks' default constructor called" << endl;
+ SetDebugLevel(1);
+ if (GetDebugLevel() > 0) cout << "AliTPCcalibTracks' default constructor called" << endl;
}
AliTPCcalibTracks::AliTPCcalibTracks(const AliTPCcalibTracks& calibTracks):
-TNamed(),
+ AliTPCcalibBase(calibTracks),
fClusterParam(0),
- fDebugStream(0),
fROC(0),
fArrayAmpRow(0),
fArrayAmp(0),
fHclusterPerPadrowRaw(0),
fClusterCutHisto(0),
fCalPadClusterPerPad(0),
- fCalPadClusterPerPadRaw(0),
- fDebugLevel(0),
- fFitterLinY1(0), //!
- fFitterLinZ1(0), //!
- fFitterLinY2(0), //!
- fFitterLinZ2(0), //!
- fFitterParY(0), //!
- fFitterParZ(0) //!
+ fCalPadClusterPerPadRaw(0)
{
//
// AliTPCcalibTracks copy constructor
//
- if (fDebugLevel > 0) cout << " ***** this is AliTPCcalibTracks' copy constructor ***** " << endl;
+ if (GetDebugLevel() > 0) cout << " ***** this is AliTPCcalibTracks' copy constructor ***** " << endl;
Bool_t dirStatus = TH1::AddDirectoryStatus();
TH1::AddDirectory(kFALSE);
fCuts = new AliTPCcalibTracksCuts(calibTracks.fCuts->GetMinClusters(), calibTracks.fCuts->GetMinRatio(),
calibTracks.fCuts->GetMax1pt(), calibTracks.fCuts->GetEdgeYXCutNoise(), calibTracks.fCuts->GetEdgeThetaCutNoise());
- fDebugLevel = calibTracks.GetLogLevel();
SetNameTitle(calibTracks.GetName(), calibTracks.GetTitle());
TH1::AddDirectory(dirStatus); // set status back to original status
// cout << "+++++ end of copy constructor +++++" << endl; // TO BE REMOVED
AliTPCcalibTracks::AliTPCcalibTracks(const Text_t *name, const Text_t *title, AliTPCClusterParam *clusterParam, AliTPCcalibTracksCuts* cuts, Int_t logLevel) :
- TNamed(name, title),
+ AliTPCcalibBase(),
fClusterParam(0),
- fDebugStream(0),
fROC(0),
fArrayAmpRow(0),
fArrayAmp(0),
fHclusterPerPadrowRaw(0),
fClusterCutHisto(0),
fCalPadClusterPerPad(0),
- fCalPadClusterPerPadRaw(0),
- fDebugLevel(0),
- fFitterLinY1(0), //!
- fFitterLinZ1(0), //!
- fFitterLinY2(0), //!
- fFitterLinZ2(0), //!
- fFitterParY(0), //!
- fFitterParZ(0) //!
+ fCalPadClusterPerPadRaw(0)
{
//
// AliTPCcalibTracks constructor
// specify 'clusterParam', (needed for TPC cluster error and shape parameterization)
// In the parameter 'cuts' the cuts are specified, that decide
// weather a track will be accepted for calibration or not.
- // log level - debug output: -1: silence, 0: default, 1: things like constructor called, 5: write fDebugStream, 6: waste your screen
+ //
+ // fDebugLevel - debug output: -1: silence, 0: default, 1: things like constructor called, 5: write fDebugStreamer, 6: waste your screen
//
// All histograms are instatiated in this constructor.
//
- if (fDebugLevel > 0) cout << " ***** this is AliTPCcalibTracks' main constructor ***** " << endl;
+ this->SetName(name);
+ this->SetTitle(title);
+
+ if (GetDebugLevel() > 0) cout << " ***** this is AliTPCcalibTracks' main constructor ***** " << endl;
G__SetCatchException(0);
fClusterParam = clusterParam;
Error("AliTPCcalibTracks","No cluster parametrization found! A valid clusterParam object is needed in the constructor. (To be found in 'TPCClusterParam.root'.)");
}
fCuts = cuts;
- fDebugLevel = logLevel;
- if (fDebugLevel > 4) fDebugStream = new TTreeSRedirector("TPCSelectorDebug.root"); // needs investigation !!!!!
+ SetDebugLevel(logLevel);
TH1::AddDirectory(kFALSE);
TProfile * prof1=0;
TH1F * his1 =0;
fHclus = new TH1I("hclus","Number of clusters per track",160, 0, 160); // valgrind 3
- fRejectedTracksHisto = new TH1I("RejectedTracksHisto", "Rejected tracks, sorted by failed cut", 10, 1, 10);
+ fRejectedTracksHisto = new TH1I("RejectedTracksHisto", "Rejected tracks, sorted by failed cut", 100, -1, 10);
fHclusterPerPadrow = new TH1I("fHclusterPerPadrow", " clusters per padRow, used for the resolution tree", 160, 0, 160);
fHclusterPerPadrowRaw = new TH1I("fHclusterPerPadrowRaw", " clusters per padRow, before cutting clusters", 160, 0, 160);
fCalPadClusterPerPad = new AliTPCCalPad("fCalPadClusterPerPad", "clusters per pad");
for (Int_t ipad = 0; ipad < 3; ipad++){
Int_t bin = GetBin(iq, ipad);
Float_t qmean = GetQ(bin);
- char name[200];
- sprintf(name,"ResolY Pad%d Qmiddle%f",ipad, qmean);
- his3D = new TH3F(name, name, 20,10,250, 20, 0,1.5, 50, -1,1);
+ char hname[200];
+ sprintf(hname,"ResolY Pad%d Qmiddle%f",ipad, qmean);
+ his3D = new TH3F(hname, hname, 20,10,250, 20, 0,1.5, 100, -1,1);
fArrayQDY->AddAt(his3D, bin);
- sprintf(name,"ResolZ Pad%d Qmiddle%f",ipad, qmean);
- his3D = new TH3F(name, name, 20,10,250, 20, 0,1.5, 50, -1,1);
+ sprintf(hname,"ResolZ Pad%d Qmiddle%f",ipad, qmean);
+ his3D = new TH3F(hname, hname, 20,10,250, 20, 0,1.5, 100, -1,1);
fArrayQDZ->AddAt(his3D, bin);
- sprintf(name,"RMSY Pad%d Qmiddle%f",ipad, qmean);
- his3D = new TH3F(name, name, 20,10,250, 20, 0,1.5, 50, 0,1);
+ sprintf(hname,"RMSY Pad%d Qmiddle%f",ipad, qmean);
+ his3D = new TH3F(hname, hname, 20,10,250, 20, 0,1.5, 100, 0,0.6);
fArrayQRMSY->AddAt(his3D, bin);
- sprintf(name,"RMSZ Pad%d Qmiddle%f",ipad, qmean);
- his3D = new TH3F(name, name, 20,10,250, 20, 0,1.5, 50, 0,1);
+ sprintf(hname,"RMSZ Pad%d Qmiddle%f",ipad, qmean);
+ his3D = new TH3F(hname, hname, 20,10,250, 20, 0,1.5, 100, 0,0.6);
fArrayQRMSZ->AddAt(his3D, bin);
}
}
}
}
- fFitterLinY1 = new TLinearFitter (2,"pol1");
- fFitterLinZ1 = new TLinearFitter (2,"pol1");
- fFitterLinY2 = new TLinearFitter (2,"pol1");
- fFitterLinZ2 = new TLinearFitter (2,"pol1");
- fFitterParY = new TLinearFitter (3,"pol2");
- fFitterParZ = new TLinearFitter (3,"pol2");
- if (fDebugLevel > 1) cout << "AliTPCcalibTracks object sucessfully constructed: " << GetName() << endl;
+ if (GetDebugLevel() > 1) cout << "AliTPCcalibTracks object sucessfully constructed: " << GetName() << endl;
cout << "end of main constructor" << endl; // TO BE REMOVED
}
// AliTPCcalibTracks destructor
//
- if (fDebugLevel > 0) cout << "AliTPCcalibTracks' destuctor called." << endl;
+ if (GetDebugLevel() > 0) cout << "AliTPCcalibTracks' destuctor called." << endl;
Int_t length = 0;
if (fArrayAmpRow) length = fArrayAmpRow->GetEntriesFast();
for (Int_t i = 0; i < length; i++){
delete fArrayChargeVsDriftlength->At(i);
}
- delete fFitterLinY1;
- delete fFitterLinZ1;
- delete fFitterLinY2;
- delete fFitterLinZ2;
- delete fFitterParY;
- delete fFitterParZ;
-
+
delete fArrayQDY;
delete fArrayQDZ;
delete fArrayQRMSY;
delete fHclusterPerPadrowRaw;
if (fCalPadClusterPerPad) delete fCalPadClusterPerPad;
if (fCalPadClusterPerPadRaw) delete fCalPadClusterPerPadRaw;
- fcalPadRegionChargeVsDriftlength->Delete();
- delete fcalPadRegionChargeVsDriftlength;
- if (fDebugLevel > 4) delete fDebugStream;
+ if(fcalPadRegionChargeVsDriftlength) {
+ fcalPadRegionChargeVsDriftlength->Delete();
+ delete fcalPadRegionChargeVsDriftlength;
+ }
}
-void AliTPCcalibTracks::AddInfo(TChain * chain, char* fileName){
- //
- // Add the neccessary information for processing to the chain
- // (cluster parametrization)
- //
- TFile clusterParamFile(fileName);
- AliTPCClusterParam *clusterParam = (AliTPCClusterParam *) clusterParamFile.Get("Param");
- chain->GetUserInfo()->AddLast((TObject*)clusterParam);
- cout << "Clusterparametrization added to the chain." << endl;
-}
-
-void AliTPCcalibTracks::Process(AliTPCseed *track, AliESDtrack */*esd*/){
+void AliTPCcalibTracks::Process(AliTPCseed *track){
//
// To be called in the selector
// first AcceptTrack is evaluated, then calls all the following analyse functions:
// FillResolutionHistoLocal(track)
// AlignUpDown(track, esd)
//
- if (fDebugLevel > 5) Info("Process","Starting to process the track...");
+ if (GetDebugLevel() > 5) Info("Process","Starting to process the track...");
Int_t accpetStatus = AcceptTrack(track);
if (accpetStatus == 0) {
FillResolutionHistoLocal(track);
//if (TMath::Abs(track->GetZ())<10.) return kFALSE;
//if (TMath::Abs(track->GetTgl())>0.03) return kFALSE;
- if (fDebugLevel > 5) Info("AcceptTrack","Track has been accepted.");
+ if (GetDebugLevel() > 20) Info("AcceptTrack","Track has been accepted.");
return 0;
}
// and to avoid redundant data
//
- if (fDebugLevel > 5) Info("FillResolutionHistoLocal"," ***** Start of FillResolutionHistoLocal *****");
+ static TLinearFitter fFitterLinY1(2,"pol1"); //
+ static TLinearFitter fFitterLinZ1(2,"pol1"); //
+ static TLinearFitter fFitterLinY2(2,"pol1"); //
+ static TLinearFitter fFitterLinZ2(2,"pol1"); //
+ static TLinearFitter fFitterParY(3,"pol2"); //
+ static TLinearFitter fFitterParZ(3,"pol2"); //
+
+ fFitterLinY1.StoreData(kFALSE);
+ fFitterLinZ1.StoreData(kFALSE);
+ fFitterLinY2.StoreData(kFALSE);
+ fFitterLinZ2.StoreData(kFALSE);
+ fFitterParY.StoreData(kFALSE);
+ fFitterParZ.StoreData(kFALSE);
+
+
+ if (GetDebugLevel() > 5) Info("FillResolutionHistoLocal"," ***** Start of FillResolutionHistoLocal *****");
const Int_t kDelta = 10; // delta rows to fit
const Float_t kMinRatio = 0.75; // minimal ratio
- const Float_t kCutChi2 = 6.; // cut chi2 - left right - kink removal
+ // const Float_t kCutChi2 = 6.; // cut chi2 - left right - kink removal
const Float_t kErrorFraction = 0.5; // use only clusters with small interpolation error - for error param
const Int_t kFirstLargePad = 127; // medium pads -> long pads
const Float_t kLargePadSize = 1.5; // factor between medium and long pads' area
const Int_t kDeltaWriteDebugStream = 5; // only for every kDeltaWriteDebugStream'th padrow debug information is calulated and written to debugstream
-// TLinearFitter fFitterLinY1 = fFitterLinY1;
-// TLinearFitter fFitterLinZ1 = ffFitterLinZ1;
-// TLinearFitter fFitterLinY2 = ffFitterLinY2;
-// TLinearFitter fFitterLinZ2 = ffFitterLinZ2;
-// TLinearFitter fFitterParY = ffFitterParY;
-// TLinearFitter fFitterParZ = ffFitterParZ;
TVectorD paramY0(2);
TVectorD paramZ0(2);
TVectorD paramY1(2);
if (sector != sectorG){
// track leaves sector before it crossed enough rows to fit / initialization
nClusters = 0;
- fFitterParY->ClearPoints();
- fFitterParZ->ClearPoints();
+ fFitterParY.ClearPoints();
+ fFitterParZ.ClearPoints();
sectorG = sector;
}
else {
nClusters++;
Double_t x = cluster0->GetX();
- fFitterParY->AddPoint(&x, cluster0->GetY(), 1);
- fFitterParZ->AddPoint(&x, cluster0->GetZ(), 1);
+ fFitterParY.AddPoint(&x, cluster0->GetY(), 1);
+ fFitterParZ.AddPoint(&x, cluster0->GetZ(), 1);
//
if ( nClusters >= kDelta + 3 ){
// if more than 13 (kDelta+3) clusters were added to the fitters
// fit the tracklet, increase trackletCounter
- fFitterParY->Eval();
- fFitterParZ->Eval();
+ fFitterParY.Eval();
+ fFitterParZ.Eval();
nTrackletsAll++;
- csigmaY += fFitterParY->GetChisquare() / (nClusters - 3.);
- csigmaZ += fFitterParZ->GetChisquare() / (nClusters - 3.);
+ csigmaY += fFitterParY.GetChisquare() / (nClusters - 3.);
+ csigmaZ += fFitterParZ.GetChisquare() / (nClusters - 3.);
nClusters = -1;
- fFitterParY->ClearPoints();
- fFitterParZ->ClearPoints();
+ fFitterParY.ClearPoints();
+ fFitterParZ.ClearPoints();
}
}
} // for (Int_t irow = 0; irow < 159; irow++)
Float_t xref = cluster0->GetX();
// Make Fit
- fFitterParY->ClearPoints();
- fFitterParZ->ClearPoints();
- fFitterLinY1->ClearPoints();
- fFitterLinZ1->ClearPoints();
- fFitterLinY2->ClearPoints();
- fFitterLinZ2->ClearPoints();
+ fFitterParY.ClearPoints();
+ fFitterParZ.ClearPoints();
+ fFitterLinY1.ClearPoints();
+ fFitterLinZ1.ClearPoints();
+ fFitterLinY2.ClearPoints();
+ fFitterLinZ2.ClearPoints();
// fit tracklet (clusters in given padrow +- kDelta padrows)
// with polynom of 2nd order and two polynoms of 1st order
nclFound++;
if (idelta < 0){
ncl0++;
- fFitterLinY1->AddPoint(&x, currentCluster->GetY(), csigmaY);
- fFitterLinZ1->AddPoint(&x, currentCluster->GetZ(), csigmaZ);
+ fFitterLinY1.AddPoint(&x, currentCluster->GetY(), csigmaY);
+ fFitterLinZ1.AddPoint(&x, currentCluster->GetZ(), csigmaZ);
}
if (idelta > 0){
ncl1++;
- fFitterLinY2->AddPoint(&x, currentCluster->GetY(), csigmaY);
- fFitterLinZ2->AddPoint(&x, currentCluster->GetZ(), csigmaZ);
+ fFitterLinY2.AddPoint(&x, currentCluster->GetY(), csigmaY);
+ fFitterLinZ2.AddPoint(&x, currentCluster->GetZ(), csigmaZ);
}
- fFitterParY->AddPoint(&x, currentCluster->GetY(), csigmaY);
- fFitterParZ->AddPoint(&x, currentCluster->GetZ(), csigmaZ);
+ fFitterParY.AddPoint(&x, currentCluster->GetY(), csigmaY);
+ fFitterParZ.AddPoint(&x, currentCluster->GetZ(), csigmaZ);
} // loop over neighbourhood for fitter filling
+
+
if (nclFound < kDelta * kMinRatio) fRejectedTracksHisto->Fill(10);
if (nclFound < kDelta * kMinRatio) fClusterCutHisto->Fill(1, irow);
if (nclFound < kDelta * kMinRatio) continue; // if not enough clusters (7.5) found in neighbourhood goto next padrow
- fFitterParY->Eval();
- fFitterParZ->Eval();
- Double_t chi2 = (fFitterParY->GetChisquare() + fFitterParZ->GetChisquare()) / (2. * nclFound - 6.);
- if (chi2 > kCutChi2) fRejectedTracksHisto->Fill(9);
- if (chi2 > kCutChi2) fClusterCutHisto->Fill(2, irow);
- if (chi2 > kCutChi2) continue; // if chi^2 is too big goto next padrow
-
+ fFitterParY.Eval();
+ fFitterParZ.Eval();
+ Double_t chi2 = (fFitterParY.GetChisquare() + fFitterParZ.GetChisquare()) / (2. * nclFound - 6.);
+ //if (chi2 > kCutChi2) fRejectedTracksHisto->Fill(9);
+ //if (chi2 > kCutChi2) fClusterCutHisto->Fill(2, irow);
+ //if (chi2 > kCutChi2) continue; // if chi^2 is too big goto next padrow
+ TTreeSRedirector *cstream = GetDebugStreamer();
+ if (cstream){
+ (*cstream)<<"Cut9"<<
+ "chi2="<<chi2<<
+ "\n";
+ }
// REMOVE KINK
// only when there are enough clusters (4) in each direction
if (ncl0 > 4){
- fFitterLinY1->Eval();
- fFitterLinZ1->Eval();
+ fFitterLinY1.Eval();
+ fFitterLinZ1.Eval();
}
if (ncl1 > 4){
- fFitterLinY2->Eval();
- fFitterLinZ2->Eval();
+ fFitterLinY2.Eval();
+ fFitterLinZ2.Eval();
}
if (ncl0 > 4 && ncl1 > 4){
- fFitterLinY1->GetCovarianceMatrix(matrixY0);
- fFitterLinY2->GetCovarianceMatrix(matrixY1);
- fFitterLinZ1->GetCovarianceMatrix(matrixZ0);
- fFitterLinZ2->GetCovarianceMatrix(matrixZ1);
- fFitterLinY2->GetParameters(paramY1);
- fFitterLinZ2->GetParameters(paramZ1);
- fFitterLinY1->GetParameters(paramY0);
- fFitterLinZ1->GetParameters(paramZ0);
+ fFitterLinY1.GetCovarianceMatrix(matrixY0);
+ fFitterLinY2.GetCovarianceMatrix(matrixY1);
+ fFitterLinZ1.GetCovarianceMatrix(matrixZ0);
+ fFitterLinZ2.GetCovarianceMatrix(matrixZ1);
+ fFitterLinY2.GetParameters(paramY1);
+ fFitterLinZ2.GetParameters(paramZ1);
+ fFitterLinY1.GetParameters(paramY0);
+ fFitterLinZ1.GetParameters(paramZ0);
paramY0 -= paramY1;
paramZ0 -= paramZ1;
matrixY0 += matrixY1;
matrixZ0 += matrixZ1;
- Double_t chi2 = 0;
+ Double_t cchi2 = 0;
TMatrixD difY(2, 1, paramY0.GetMatrixArray());
TMatrixD difYT(1, 2, paramY0.GetMatrixArray());
matrixY0.Invert();
TMatrixD mulY(matrixY0, TMatrixD::kMult, difY);
TMatrixD chi2Y(difYT, TMatrixD::kMult, mulY);
- chi2 += chi2Y(0, 0);
+ cchi2 += chi2Y(0, 0);
TMatrixD difZ(2, 1, paramZ0.GetMatrixArray());
TMatrixD difZT(1, 2, paramZ0.GetMatrixArray());
matrixZ0.Invert();
TMatrixD mulZ(matrixZ0, TMatrixD::kMult, difZ);
TMatrixD chi2Z(difZT, TMatrixD::kMult, mulZ);
- chi2 += chi2Z(0, 0);
+ cchi2 += chi2Z(0, 0);
- // REMOVE KINK
- if (chi2 * 0.25 > kCutChi2) fRejectedTracksHisto->Fill(8);
- if (chi2 * 0.25 > kCutChi2) fClusterCutHisto->Fill(3, irow);
- if (chi2 * 0.25 > kCutChi2) continue; // if chi2 is too big goto next padrow
+ // REMOVE KINK - TO be fixed - proper chi2 calculation for curved track to be implemented
+ //if (chi2 * 0.25 > kCutChi2) fRejectedTracksHisto->Fill(8);
+ //if (chi2 * 0.25 > kCutChi2) fClusterCutHisto->Fill(3, irow);
+ //if (chi2 * 0.25 > kCutChi2) continue; // if chi2 is too big goto next padrow
// fit tracklet with polynom of 2nd order and two polynoms of 1st order
// take both polynoms of 1st order, calculate difference of their parameters
// add covariance matrixes and calculate chi2 of this difference
// if this chi2 is bigger than a given threshold, assume that the current cluster is
// a kink an goto next padrow
+
+ if (cstream){
+ (*cstream)<<"Cut8"<<
+ "chi2="<<cchi2<<
+ "\n";
+ }
}
// current padrow has no kink
// get fit parameters from pol2 fit:
Double_t paramY[4], paramZ[4];
- paramY[0] = fFitterParY->GetParameter(0);
- paramY[1] = fFitterParY->GetParameter(1);
- paramY[2] = fFitterParY->GetParameter(2);
- paramZ[0] = fFitterParZ->GetParameter(0);
- paramZ[1] = fFitterParZ->GetParameter(1);
- paramZ[2] = fFitterParZ->GetParameter(2);
+ paramY[0] = fFitterParY.GetParameter(0);
+ paramY[1] = fFitterParY.GetParameter(1);
+ paramY[2] = fFitterParY.GetParameter(2);
+ paramZ[0] = fFitterParZ.GetParameter(0);
+ paramZ[1] = fFitterParZ.GetParameter(1);
+ paramZ[2] = fFitterParZ.GetParameter(2);
Double_t tracky = paramY[0];
Double_t trackz = paramZ[0];
// Fill resolution histograms
Bool_t useForResol = kTRUE;
- if (fFitterParY->GetParError(0) > kErrorFraction * csigmaY) useForResol = kFALSE;
-
+ if (fFitterParY.GetParError(0) > kErrorFraction * csigmaY) useForResol = kFALSE;
+
+ if (cstream){
+ Float_t zdrift = 250 - TMath::Abs(cluster0->GetZ());
+ Float_t sy = cluster0->GetSigmaY2();
+ Float_t sz = cluster0->GetSigmaZ2();
+ (*cstream)<<"Resol0"<<
+ "run="<<fRun<< // run number
+ "event="<<fEvent<< // event number
+ "time="<<fTime<< // time stamp of event
+ "trigger="<<fTrigger<< // trigger
+ "mag="<<fMagF<< // magnetic field
+ "padSize="<<padSize<<
+ "angley="<<angley<<
+ "anglez="<<anglez<<
+ "zdr="<<zdrift<<
+ "dy="<<deltay<<
+ "dz="<<deltaz<<
+ "sy="<<sy<<
+ "sz="<<sz<<
+ "\n";
+ }
+
if (useForResol){
fDeltaY->Fill(deltay);
fDeltaZ->Fill(deltaz);
//=============================================================================================
if (useForResol && nclFound > 2 * kMinRatio * kDelta
- && irow % kDeltaWriteDebugStream == 0 && fDebugLevel > 4){
- if (fDebugLevel > 5) Info("FillResolutionHistoLocal","Filling 'TPCSelectorDebug.root', irow = %i", irow);
+ && irow % kDeltaWriteDebugStream == 0 && GetDebugLevel() > 4){
+ if (GetDebugLevel() > 20) Info("FillResolutionHistoLocal","Filling 'TPCSelectorDebug.root', irow = %i", irow);
FillResolutionHistoLocalDebugPart(track, cluster0, irow, angley, anglez, nclFound, kDelta);
} // if (useForResol && nclFound > 2 * kMinRatio * kDelta)
//
// - debug part of FillResolutionHistoLocal -
// called only for every kDeltaWriteDebugStream'th padrow, to avoid to much redundant data
- // called only for fDebugLevel > 4
+ // called only for GetStreamLevel() > 4
// fill resolution trees
//
if (cluster0->GetRow() > 63) padSize = 2; // long pads
}
- static TLinearFitter fitY0(3, "pol2");
- static TLinearFitter fitZ0(3, "pol2");
- static TLinearFitter fitY2(5, "hyp4");
- static TLinearFitter fitZ2(5, "hyp4");
- static TLinearFitter fitY2Q(5, "hyp4");
- static TLinearFitter fitZ2Q(5, "hyp4");
- static TLinearFitter fitY2S(5, "hyp4");
- static TLinearFitter fitZ2S(5, "hyp4");
- fitY0.ClearPoints();
- fitZ0.ClearPoints();
- fitY2.ClearPoints();
- fitZ2.ClearPoints();
- fitY2Q.ClearPoints();
- fitZ2Q.ClearPoints();
- fitY2S.ClearPoints();
- fitZ2S.ClearPoints();
-
- for (Int_t idelta = -kDelta; idelta <= kDelta; idelta++){
- // loop over irow +- kDelta rows (neighboured rows)
- //
- //
- if (idelta == 0) continue;
- if (idelta + irow < 0 || idelta + irow > 159) continue; // don't go out of ROC
- AliTPCclusterMI * cluster = track->GetClusterPointer(irow + idelta);
- if (!cluster) continue;
- if (cluster->GetType() < 0) continue;
- if (cluster->GetDetector() != sector) continue;
- Double_t x = cluster->GetX() - xref;
- Double_t sigmaY0 = fClusterParam->GetError0Par( 0, padSize, (250.0 - TMath::Abs(cluster->GetZ())), TMath::Abs(angley) );
- Double_t sigmaZ0 = fClusterParam->GetError0Par( 1, padSize, (250.0 - TMath::Abs(cluster->GetZ())), TMath::Abs(anglez) );
- //
- Double_t sigmaYQ = fClusterParam->GetErrorQPar( 0, padSize, (250.0 - TMath::Abs(cluster->GetZ())), TMath::Abs(angley), TMath::Abs(cluster->GetMax()) );
- Double_t sigmaZQ = fClusterParam->GetErrorQPar( 1, padSize, (250.0 - TMath::Abs(cluster->GetZ())), TMath::Abs(anglez), TMath::Abs(cluster->GetMax()) );
- Double_t sigmaYS = fClusterParam->GetErrorQParScaled( 0, padSize, (250.0 - TMath::Abs(cluster->GetZ())),
+ static TLinearFitter fitY0(3, "pol2");
+ static TLinearFitter fitZ0(3, "pol2");
+ static TLinearFitter fitY2(5, "hyp4");
+ static TLinearFitter fitZ2(5, "hyp4");
+ static TLinearFitter fitY2Q(5, "hyp4");
+ static TLinearFitter fitZ2Q(5, "hyp4");
+ static TLinearFitter fitY2S(5, "hyp4");
+ static TLinearFitter fitZ2S(5, "hyp4");
+ fitY0.ClearPoints();
+ fitZ0.ClearPoints();
+ fitY2.ClearPoints();
+ fitZ2.ClearPoints();
+ fitY2Q.ClearPoints();
+ fitZ2Q.ClearPoints();
+ fitY2S.ClearPoints();
+ fitZ2S.ClearPoints();
+
+ for (Int_t idelta = -kDelta; idelta <= kDelta; idelta++){
+ // loop over irow +- kDelta rows (neighboured rows)
+ //
+ //
+ if (idelta == 0) continue;
+ if (idelta + irow < 0 || idelta + irow > 159) continue; // don't go out of ROC
+ AliTPCclusterMI * cluster = track->GetClusterPointer(irow + idelta);
+ if (!cluster) continue;
+ if (cluster->GetType() < 0) continue;
+ if (cluster->GetDetector() != sector) continue;
+ Double_t x = cluster->GetX() - xref;
+ Double_t sigmaY0 = fClusterParam->GetError0Par( 0, padSize, (250.0 - TMath::Abs(cluster->GetZ())), TMath::Abs(angley) );
+ Double_t sigmaZ0 = fClusterParam->GetError0Par( 1, padSize, (250.0 - TMath::Abs(cluster->GetZ())), TMath::Abs(anglez) );
+ //
+ Double_t sigmaYQ = fClusterParam->GetErrorQPar( 0, padSize, (250.0 - TMath::Abs(cluster->GetZ())), TMath::Abs(angley), TMath::Abs(cluster->GetMax()) );
+ Double_t sigmaZQ = fClusterParam->GetErrorQPar( 1, padSize, (250.0 - TMath::Abs(cluster->GetZ())), TMath::Abs(anglez), TMath::Abs(cluster->GetMax()) );
+ Double_t sigmaYS = fClusterParam->GetErrorQParScaled( 0, padSize, (250.0 - TMath::Abs(cluster->GetZ())),
TMath::Abs(angley), TMath::Abs(cluster->GetMax()) );
- Double_t sigmaZS = fClusterParam->GetErrorQParScaled( 1, padSize, (250.0 - TMath::Abs(cluster->GetZ())),
+ Double_t sigmaZS = fClusterParam->GetErrorQParScaled( 1, padSize, (250.0 - TMath::Abs(cluster->GetZ())),
TMath::Abs(anglez), TMath::Abs(cluster->GetMax()) );
- Float_t rmsYFactor = fClusterParam->GetShapeFactor( 0, padSize,(250.0 - TMath::Abs(cluster->GetZ())),
- TMath::Abs(anglez), TMath::Abs(cluster->GetMax()),
- TMath::Sqrt(cluster0->GetSigmaY2()), 0 );
- Float_t rmsZFactor = fClusterParam->GetShapeFactor(0, padSize,(250.0 - TMath::Abs(cluster->GetZ())),
- TMath::Abs(anglez), TMath::Abs(cluster->GetMax()),
- TMath::Sqrt(cluster0->GetSigmaZ2()),0 );
- sigmaYS = TMath::Sqrt(sigmaYS * sigmaYS + rmsYFactor * rmsYFactor / 12.);
- sigmaZS = TMath::Sqrt(sigmaZS * sigmaZS + rmsZFactor * rmsZFactor / 12. + rmsYFactor * rmsYFactor / 24.);
- //
- if (kDelta != 0){
- fitY0.AddPoint(&x, cluster->GetY(), sigmaY0);
- fitZ0.AddPoint(&x, cluster->GetZ(), sigmaZ0);
- }
- Double_t xxx[4];
- xxx[0] = ( (idelta+irow) % 2 == 0 ) ? 1 : 0;
- xxx[1] = x;
- xxx[2] = ( (idelta+irow) % 2 == 0 ) ? x : 0;
- xxx[3] = x * x;
- fitY2.AddPoint(xxx, cluster->GetY(), sigmaY0);
- fitY2Q.AddPoint(xxx, cluster->GetY(), sigmaYQ);
- fitY2S.AddPoint(xxx, cluster->GetY(), sigmaYS);
- fitZ2.AddPoint(xxx, cluster->GetZ(), sigmaZ0);
- fitZ2Q.AddPoint(xxx, cluster->GetZ(), sigmaZQ);
- fitZ2S.AddPoint(xxx, cluster->GetZ(), sigmaZS);
- //
- } // neigbouhood-loop
- //
- fitY0.Eval();
- fitZ0.Eval();
- fitY2.Eval();
- fitZ2.Eval();
- fitY2Q.Eval();
- fitZ2Q.Eval();
- fitY2S.Eval();
- fitZ2S.Eval();
- Float_t chi2Y0 = fitY0.GetChisquare() / (nclFound-3.);
- Float_t chi2Z0 = fitZ0.GetChisquare() / (nclFound-3.);
- Float_t chi2Y2 = fitY2.GetChisquare() / (nclFound-5.);
- Float_t chi2Z2 = fitZ2.GetChisquare() / (nclFound-5.);
- Float_t chi2Y2Q = fitY2Q.GetChisquare() / (nclFound-5.);
- Float_t chi2Z2Q = fitZ2Q.GetChisquare() / (nclFound-5.);
- Float_t chi2Y2S = fitY2S.GetChisquare() / (nclFound-5.);
- Float_t chi2Z2S = fitZ2S.GetChisquare() / (nclFound-5.);
- //
- static TVectorD parY0(3);
- static TMatrixD matY0(3, 3);
- static TVectorD parZ0(3);
- static TMatrixD matZ0(3, 3);
- fitY0.GetParameters(parY0);
- fitY0.GetCovarianceMatrix(matY0);
- fitZ0.GetParameters(parZ0);
- fitZ0.GetCovarianceMatrix(matZ0);
- //
- static TVectorD parY2(5);
- static TMatrixD matY2(5,5);
- static TVectorD parZ2(5);
- static TMatrixD matZ2(5,5);
- fitY2.GetParameters(parY2);
- fitY2.GetCovarianceMatrix(matY2);
- fitZ2.GetParameters(parZ2);
- fitZ2.GetCovarianceMatrix(matZ2);
- //
- static TVectorD parY2Q(5);
- static TMatrixD matY2Q(5,5);
- static TVectorD parZ2Q(5);
- static TMatrixD matZ2Q(5,5);
- fitY2Q.GetParameters(parY2Q);
- fitY2Q.GetCovarianceMatrix(matY2Q);
- fitZ2Q.GetParameters(parZ2Q);
- fitZ2Q.GetCovarianceMatrix(matZ2Q);
- static TVectorD parY2S(5);
- static TMatrixD matY2S(5,5);
- static TVectorD parZ2S(5);
- static TMatrixD matZ2S(5,5);
- fitY2S.GetParameters(parY2S);
- fitY2S.GetCovarianceMatrix(matY2S);
- fitZ2S.GetParameters(parZ2S);
- fitZ2S.GetCovarianceMatrix(matZ2S);
- Float_t sigmaY0 = TMath::Sqrt(matY0(0,0));
- Float_t sigmaZ0 = TMath::Sqrt(matZ0(0,0));
- Float_t sigmaDY0 = TMath::Sqrt(matY0(1,1));
- Float_t sigmaDZ0 = TMath::Sqrt(matZ0(1,1));
- Float_t sigmaY2 = TMath::Sqrt(matY2(1,1));
- Float_t sigmaZ2 = TMath::Sqrt(matZ2(1,1));
- Float_t sigmaDY2 = TMath::Sqrt(matY2(3,3));
- Float_t sigmaDZ2 = TMath::Sqrt(matZ2(3,3));
- Float_t sigmaY2Q = TMath::Sqrt(matY2Q(1,1));
- Float_t sigmaZ2Q = TMath::Sqrt(matZ2Q(1,1));
- Float_t sigmaDY2Q = TMath::Sqrt(matY2Q(3,3));
- Float_t sigmaDZ2Q = TMath::Sqrt(matZ2Q(3,3));
- Float_t sigmaY2S = TMath::Sqrt(matY2S(1,1));
- Float_t sigmaZ2S = TMath::Sqrt(matZ2S(1,1));
- Float_t sigmaDY2S = TMath::Sqrt(matY2S(3,3));
- Float_t sigmaDZ2S = TMath::Sqrt(matZ2S(3,3));
-
- // Error parameters
- Float_t csigmaY0 = fClusterParam->GetError0Par(0,padSize,(250.0-TMath::Abs(cluster0->GetZ())),TMath::Abs(angley));
- Float_t csigmaZ0 = fClusterParam->GetError0Par(1,padSize,(250.0-TMath::Abs(cluster0->GetZ())),TMath::Abs(anglez));
+ Float_t rmsYFactor = fClusterParam->GetShapeFactor( 0, padSize,(250.0 - TMath::Abs(cluster->GetZ())),
+ TMath::Abs(anglez), TMath::Abs(cluster->GetMax()),
+ TMath::Sqrt(cluster0->GetSigmaY2()), 0 );
+ Float_t rmsZFactor = fClusterParam->GetShapeFactor(0, padSize,(250.0 - TMath::Abs(cluster->GetZ())),
+ TMath::Abs(anglez), TMath::Abs(cluster->GetMax()),
+ TMath::Sqrt(cluster0->GetSigmaZ2()),0 );
+ sigmaYS = TMath::Sqrt(sigmaYS * sigmaYS + rmsYFactor * rmsYFactor / 12.);
+ sigmaZS = TMath::Sqrt(sigmaZS * sigmaZS + rmsZFactor * rmsZFactor / 12. + rmsYFactor * rmsYFactor / 24.);
+ //
+ if (kDelta != 0){
+ fitY0.AddPoint(&x, cluster->GetY(), sigmaY0);
+ fitZ0.AddPoint(&x, cluster->GetZ(), sigmaZ0);
+ }
+ Double_t xxx[4];
+ xxx[0] = ( (idelta+irow) % 2 == 0 ) ? 1 : 0;
+ xxx[1] = x;
+ xxx[2] = ( (idelta+irow) % 2 == 0 ) ? x : 0;
+ xxx[3] = x * x;
+ fitY2.AddPoint(xxx, cluster->GetY(), sigmaY0);
+ fitY2Q.AddPoint(xxx, cluster->GetY(), sigmaYQ);
+ fitY2S.AddPoint(xxx, cluster->GetY(), sigmaYS);
+ fitZ2.AddPoint(xxx, cluster->GetZ(), sigmaZ0);
+ fitZ2Q.AddPoint(xxx, cluster->GetZ(), sigmaZQ);
+ fitZ2S.AddPoint(xxx, cluster->GetZ(), sigmaZS);
+ //
+ } // neigbouhood-loop
//
- Float_t csigmaYQ = fClusterParam->GetErrorQPar(0,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
- TMath::Abs(angley), TMath::Abs(cluster0->GetMax()));
- Float_t csigmaZQ = fClusterParam->GetErrorQPar(1,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
- TMath::Abs(anglez),TMath::Abs(cluster0->GetMax()));
- Float_t csigmaYS = fClusterParam->GetErrorQParScaled(0,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
- TMath::Abs(angley), TMath::Abs(cluster0->GetMax()));
- Float_t csigmaZS = fClusterParam->GetErrorQParScaled(1,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
- TMath::Abs(anglez),TMath::Abs(cluster0->GetMax()));
- ///////////////////////////////////////////////////////////////////////////////
-// //
-// Class to analyse tracks for calibration //
-// to be used as a component in AliTPCSelectorTracks //
-// In the constructor you have to specify name and title //
-// to get the Object out of a file. //
-// The parameter 'clusterParam', a AliTPCClusterParam object //
-// (needed for TPC cluster error and shape parameterization) //
-// Normally you get this object out of the file 'TPCClusterParam.root' //
-// In the parameter 'cuts' the cuts are specified, that decide //
-// weather a track will be accepted for calibration or not. //
-// //
-// //
-///////////////////////////////////////////////////////////////////////////////
-
- // RMS parameters
- Float_t meanRMSY = 0;
- Float_t meanRMSZ = 0;
- Int_t nclRMS = 0;
- for (Int_t idelta = -2; idelta <= 2; idelta++){
- // loop over neighbourhood
- if (idelta+irow < 0 || idelta+irow > 159) continue;
-// if (idelta+irow>159) continue;
- AliTPCclusterMI * cluster = track->GetClusterPointer(irow+idelta);
- if (!cluster) continue;
- meanRMSY += TMath::Sqrt(cluster->GetSigmaY2());
- meanRMSZ += TMath::Sqrt(cluster->GetSigmaZ2());
- nclRMS++;
- }
- meanRMSY /= nclRMS;
- meanRMSZ /= nclRMS;
-
- Float_t rmsY = TMath::Sqrt(cluster0->GetSigmaY2());
- Float_t rmsZ = TMath::Sqrt(cluster0->GetSigmaZ2());
- Float_t rmsYT = fClusterParam->GetRMSQ(0,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
- TMath::Abs(angley), TMath::Abs(cluster0->GetMax()));
- Float_t rmsZT = fClusterParam->GetRMSQ(1,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
- TMath::Abs(anglez), TMath::Abs(cluster0->GetMax()));
- Float_t rmsYT0 = fClusterParam->GetRMS0(0,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
- TMath::Abs(angley));
- Float_t rmsZT0 = fClusterParam->GetRMS0(1,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
+ fitY0.Eval();
+ fitZ0.Eval();
+ fitY2.Eval();
+ fitZ2.Eval();
+ fitY2Q.Eval();
+ fitZ2Q.Eval();
+ fitY2S.Eval();
+ fitZ2S.Eval();
+ Float_t chi2Y0 = fitY0.GetChisquare() / (nclFound-3.);
+ Float_t chi2Z0 = fitZ0.GetChisquare() / (nclFound-3.);
+ Float_t chi2Y2 = fitY2.GetChisquare() / (nclFound-5.);
+ Float_t chi2Z2 = fitZ2.GetChisquare() / (nclFound-5.);
+ Float_t chi2Y2Q = fitY2Q.GetChisquare() / (nclFound-5.);
+ Float_t chi2Z2Q = fitZ2Q.GetChisquare() / (nclFound-5.);
+ Float_t chi2Y2S = fitY2S.GetChisquare() / (nclFound-5.);
+ Float_t chi2Z2S = fitZ2S.GetChisquare() / (nclFound-5.);
+ //
+ static TVectorD parY0(3);
+ static TMatrixD matY0(3, 3);
+ static TVectorD parZ0(3);
+ static TMatrixD matZ0(3, 3);
+ fitY0.GetParameters(parY0);
+ fitY0.GetCovarianceMatrix(matY0);
+ fitZ0.GetParameters(parZ0);
+ fitZ0.GetCovarianceMatrix(matZ0);
+ //
+ static TVectorD parY2(5);
+ static TMatrixD matY2(5,5);
+ static TVectorD parZ2(5);
+ static TMatrixD matZ2(5,5);
+ fitY2.GetParameters(parY2);
+ fitY2.GetCovarianceMatrix(matY2);
+ fitZ2.GetParameters(parZ2);
+ fitZ2.GetCovarianceMatrix(matZ2);
+ //
+ static TVectorD parY2Q(5);
+ static TMatrixD matY2Q(5,5);
+ static TVectorD parZ2Q(5);
+ static TMatrixD matZ2Q(5,5);
+ fitY2Q.GetParameters(parY2Q);
+ fitY2Q.GetCovarianceMatrix(matY2Q);
+ fitZ2Q.GetParameters(parZ2Q);
+ fitZ2Q.GetCovarianceMatrix(matZ2Q);
+ static TVectorD parY2S(5);
+ static TMatrixD matY2S(5,5);
+ static TVectorD parZ2S(5);
+ static TMatrixD matZ2S(5,5);
+ fitY2S.GetParameters(parY2S);
+ fitY2S.GetCovarianceMatrix(matY2S);
+ fitZ2S.GetParameters(parZ2S);
+ fitZ2S.GetCovarianceMatrix(matZ2S);
+ Float_t sigmaY0 = TMath::Sqrt(matY0(0,0));
+ Float_t sigmaZ0 = TMath::Sqrt(matZ0(0,0));
+ Float_t sigmaDY0 = TMath::Sqrt(matY0(1,1));
+ Float_t sigmaDZ0 = TMath::Sqrt(matZ0(1,1));
+ Float_t sigmaY2 = TMath::Sqrt(matY2(1,1));
+ Float_t sigmaZ2 = TMath::Sqrt(matZ2(1,1));
+ Float_t sigmaDY2 = TMath::Sqrt(matY2(3,3));
+ Float_t sigmaDZ2 = TMath::Sqrt(matZ2(3,3));
+ Float_t sigmaY2Q = TMath::Sqrt(matY2Q(1,1));
+ Float_t sigmaZ2Q = TMath::Sqrt(matZ2Q(1,1));
+ Float_t sigmaDY2Q = TMath::Sqrt(matY2Q(3,3));
+ Float_t sigmaDZ2Q = TMath::Sqrt(matZ2Q(3,3));
+ Float_t sigmaY2S = TMath::Sqrt(matY2S(1,1));
+ Float_t sigmaZ2S = TMath::Sqrt(matZ2S(1,1));
+ Float_t sigmaDY2S = TMath::Sqrt(matY2S(3,3));
+ Float_t sigmaDZ2S = TMath::Sqrt(matZ2S(3,3));
+
+ // Error parameters
+ Float_t csigmaY0 = fClusterParam->GetError0Par(0,padSize,(250.0-TMath::Abs(cluster0->GetZ())),TMath::Abs(angley));
+ Float_t csigmaZ0 = fClusterParam->GetError0Par(1,padSize,(250.0-TMath::Abs(cluster0->GetZ())),TMath::Abs(anglez));
+ //
+ Float_t csigmaYQ = fClusterParam->GetErrorQPar(0,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
+ TMath::Abs(angley), TMath::Abs(cluster0->GetMax()));
+ Float_t csigmaZQ = fClusterParam->GetErrorQPar(1,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
+ TMath::Abs(anglez),TMath::Abs(cluster0->GetMax()));
+ Float_t csigmaYS = fClusterParam->GetErrorQParScaled(0,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
+ TMath::Abs(angley), TMath::Abs(cluster0->GetMax()));
+ Float_t csigmaZS = fClusterParam->GetErrorQParScaled(1,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
+ TMath::Abs(anglez),TMath::Abs(cluster0->GetMax()));
+
+
+ // RMS parameters
+ Float_t meanRMSY = 0;
+ Float_t meanRMSZ = 0;
+ Int_t nclRMS = 0;
+ for (Int_t idelta = -2; idelta <= 2; idelta++){
+ // loop over neighbourhood
+ if (idelta+irow < 0 || idelta+irow > 159) continue;
+ // if (idelta+irow>159) continue;
+ AliTPCclusterMI * cluster = track->GetClusterPointer(irow+idelta);
+ if (!cluster) continue;
+ meanRMSY += TMath::Sqrt(cluster->GetSigmaY2());
+ meanRMSZ += TMath::Sqrt(cluster->GetSigmaZ2());
+ nclRMS++;
+ }
+ meanRMSY /= nclRMS;
+ meanRMSZ /= nclRMS;
+
+ Float_t rmsY = TMath::Sqrt(cluster0->GetSigmaY2());
+ Float_t rmsZ = TMath::Sqrt(cluster0->GetSigmaZ2());
+ Float_t rmsYT = fClusterParam->GetRMSQ(0,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
+ TMath::Abs(angley), TMath::Abs(cluster0->GetMax()));
+ Float_t rmsZT = fClusterParam->GetRMSQ(1,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
+ TMath::Abs(anglez), TMath::Abs(cluster0->GetMax()));
+ Float_t rmsYT0 = fClusterParam->GetRMS0(0,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
+ TMath::Abs(angley));
+ Float_t rmsZT0 = fClusterParam->GetRMS0(1,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
TMath::Abs(anglez));
- Float_t rmsYSigma = fClusterParam->GetRMSSigma(0,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
- TMath::Abs(anglez), TMath::Abs(cluster0->GetMax()));
- Float_t rmsZSigma = fClusterParam->GetRMSSigma(0,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
- TMath::Abs(anglez), TMath::Abs(cluster0->GetMax()));
- Float_t rmsYFactor = fClusterParam->GetShapeFactor(0,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
- TMath::Abs(anglez), TMath::Abs(cluster0->GetMax()),
- rmsY,meanRMSY);
- Float_t rmsZFactor = fClusterParam->GetShapeFactor(0,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
- TMath::Abs(anglez), TMath::Abs(cluster0->GetMax()),
- rmsZ,meanRMSZ);
-
- // cluster debug
- (*fDebugStream)<<"ResolCl"<< // valgrind 3 40,000 bytes in 1 blocks are possibly
- "Sector="<<sector<<
- "Cl.="<<cluster0<<
- "CSigmaY0="<<csigmaY0<< // cluster errorY
- "CSigmaYQ="<<csigmaYQ<< // cluster errorY - q scaled
- "CSigmaYS="<<csigmaYS<< // cluster errorY - q scaled
- "CSigmaZ0="<<csigmaZ0<< //
- "CSigmaZQ="<<csigmaZQ<<
- "CSigmaZS="<<csigmaZS<<
- "shapeYF="<<rmsYFactor<<
- "shapeZF="<<rmsZFactor<<
- "rmsY="<<rmsY<<
- "rmsZ="<<rmsZ<<
- "rmsYM="<<meanRMSY<<
- "rmsZM="<<meanRMSZ<<
- "rmsYT="<<rmsYT<<
- "rmsZT="<<rmsZT<<
- "rmsYT0="<<rmsYT0<<
- "rmsZT0="<<rmsZT0<<
- "rmsYS="<<rmsYSigma<<
- "rmsZS="<<rmsZSigma<<
- "padSize="<<padSize<<
- "Ncl="<<nclFound<<
- "PY0.="<<&parY0<<
- "PZ0.="<<&parZ0<<
- "SigmaY0="<<sigmaY0<<
- "SigmaZ0="<<sigmaZ0<<
- "angley="<<angley<<
- "anglez="<<anglez<<
- "\n";
-
-// tracklet dubug
- (*fDebugStream)<<"ResolTr"<<
- "padSize="<<padSize<<
- "IPad="<<padSize<<
- "Sector="<<sector<<
- "Ncl="<<nclFound<<
- "chi2Y0="<<chi2Y0<<
- "chi2Z0="<<chi2Z0<<
- "chi2Y2="<<chi2Y2<<
- "chi2Z2="<<chi2Z2<<
- "chi2Y2Q="<<chi2Y2Q<<
- "chi2Z2Q="<<chi2Z2Q<<
- "chi2Y2S="<<chi2Y2S<<
- "chi2Z2S="<<chi2Z2S<<
- "PY0.="<<&parY0<<
- "PZ0.="<<&parZ0<<
- "PY2.="<<&parY2<<
- "PZ2.="<<&parZ2<<
- "PY2Q.="<<&parY2Q<<
- "PZ2Q.="<<&parZ2Q<<
- "PY2S.="<<&parY2S<<
- "PZ2S.="<<&parZ2S<<
- "SigmaY0="<<sigmaY0<<
- "SigmaZ0="<<sigmaZ0<<
- "SigmaDY0="<<sigmaDY0<<
- "SigmaDZ0="<<sigmaDZ0<<
- "SigmaY2="<<sigmaY2<<
- "SigmaZ2="<<sigmaZ2<<
- "SigmaDY2="<<sigmaDY2<<
- "SigmaDZ2="<<sigmaDZ2<<
- "SigmaY2Q="<<sigmaY2Q<<
- "SigmaZ2Q="<<sigmaZ2Q<<
- "SigmaDY2Q="<<sigmaDY2Q<<
- "SigmaDZ2Q="<<sigmaDZ2Q<<
- "SigmaY2S="<<sigmaY2S<<
- "SigmaZ2S="<<sigmaZ2S<<
- "SigmaDY2S="<<sigmaDY2S<<
- "SigmaDZ2S="<<sigmaDZ2S<<
- "angley="<<angley<<
- "anglez="<<anglez<<
- "\n";
-
-
+ Float_t rmsYSigma = fClusterParam->GetRMSSigma(0,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
+ TMath::Abs(anglez), TMath::Abs(cluster0->GetMax()));
+ Float_t rmsZSigma = fClusterParam->GetRMSSigma(0,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
+ TMath::Abs(anglez), TMath::Abs(cluster0->GetMax()));
+ Float_t rmsYFactor = fClusterParam->GetShapeFactor(0,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
+ TMath::Abs(anglez), TMath::Abs(cluster0->GetMax()),
+ rmsY,meanRMSY);
+ Float_t rmsZFactor = fClusterParam->GetShapeFactor(0,padSize,(250.0-TMath::Abs(cluster0->GetZ())),
+ TMath::Abs(anglez), TMath::Abs(cluster0->GetMax()),
+ rmsZ,meanRMSZ);
+ //
+ // cluster debug
+ TTreeSRedirector *cstream = GetDebugStreamer();
+ if (cstream){
+ (*cstream)<<"ResolCl"<< // valgrind 3 40,000 bytes in 1 blocks are possibly
+ "run="<<fRun<< // run number
+ "event="<<fEvent<< // event number
+ "time="<<fTime<< // time stamp of event
+ "trigger="<<fTrigger<< // trigger
+ "mag="<<fMagF<< // magnetic field
+ "Sector="<<sector<<
+ "Cl.="<<cluster0<<
+ "CSigmaY0="<<csigmaY0<< // cluster errorY
+ "CSigmaYQ="<<csigmaYQ<< // cluster errorY - q scaled
+ "CSigmaYS="<<csigmaYS<< // cluster errorY - q scaled
+ "CSigmaZ0="<<csigmaZ0<< //
+ "CSigmaZQ="<<csigmaZQ<<
+ "CSigmaZS="<<csigmaZS<<
+ "shapeYF="<<rmsYFactor<<
+ "shapeZF="<<rmsZFactor<<
+ "rmsY="<<rmsY<<
+ "rmsZ="<<rmsZ<<
+ "rmsYM="<<meanRMSY<<
+ "rmsZM="<<meanRMSZ<<
+ "rmsYT="<<rmsYT<<
+ "rmsZT="<<rmsZT<<
+ "rmsYT0="<<rmsYT0<<
+ "rmsZT0="<<rmsZT0<<
+ "rmsYS="<<rmsYSigma<<
+ "rmsZS="<<rmsZSigma<<
+ "padSize="<<padSize<<
+ "Ncl="<<nclFound<<
+ "PY0.="<<&parY0<<
+ "PZ0.="<<&parZ0<<
+ "SigmaY0="<<sigmaY0<<
+ "SigmaZ0="<<sigmaZ0<<
+ "angley="<<angley<<
+ "anglez="<<anglez<<
+ "\n";
+
+ // tracklet dubug
+ (*cstream)<<"ResolTr"<<
+ "run="<<fRun<< // run number
+ "event="<<fEvent<< // event number
+ "time="<<fTime<< // time stamp of event
+ "trigger="<<fTrigger<< // trigger
+ "mag="<<fMagF<< // magnetic field
+ "padSize="<<padSize<<
+ "IPad="<<padSize<<
+ "Sector="<<sector<<
+ "Ncl="<<nclFound<<
+ "chi2Y0="<<chi2Y0<<
+ "chi2Z0="<<chi2Z0<<
+ "chi2Y2="<<chi2Y2<<
+ "chi2Z2="<<chi2Z2<<
+ "chi2Y2Q="<<chi2Y2Q<<
+ "chi2Z2Q="<<chi2Z2Q<<
+ "chi2Y2S="<<chi2Y2S<<
+ "chi2Z2S="<<chi2Z2S<<
+ "PY0.="<<&parY0<<
+ "PZ0.="<<&parZ0<<
+ "PY2.="<<&parY2<<
+ "PZ2.="<<&parZ2<<
+ "PY2Q.="<<&parY2Q<<
+ "PZ2Q.="<<&parZ2Q<<
+ "PY2S.="<<&parY2S<<
+ "PZ2S.="<<&parZ2S<<
+ "SigmaY0="<<sigmaY0<<
+ "SigmaZ0="<<sigmaZ0<<
+ "SigmaDY0="<<sigmaDY0<<
+ "SigmaDZ0="<<sigmaDZ0<<
+ "SigmaY2="<<sigmaY2<<
+ "SigmaZ2="<<sigmaZ2<<
+ "SigmaDY2="<<sigmaDY2<<
+ "SigmaDZ2="<<sigmaDZ2<<
+ "SigmaY2Q="<<sigmaY2Q<<
+ "SigmaZ2Q="<<sigmaZ2Q<<
+ "SigmaDY2Q="<<sigmaDY2Q<<
+ "SigmaDZ2Q="<<sigmaDZ2Q<<
+ "SigmaY2S="<<sigmaY2S<<
+ "SigmaZ2S="<<sigmaZ2S<<
+ "SigmaDY2S="<<sigmaDY2S<<
+ "SigmaDZ2S="<<sigmaDZ2S<<
+ "angley="<<angley<<
+ "anglez="<<anglez<<
+ "\n";
+ }
}
// will draws some exemplaric pictures
//
- if (fDebugLevel > 6) Info("Draw", "Drawing an exemplaric picture.");
+ if (GetDebugLevel() > 6) Info("Draw", "Drawing an exemplaric picture.");
SetStyle();
Double_t min = 0;
Double_t max = 0;
}
-void AliTPCcalibTracks::MakeReport(Int_t stat, char* pathName){
+void AliTPCcalibTracks::MakeReport(Int_t stat, const char* pathName){
//
// all functions are called, that produce pictures
// the histograms are written to the directory 'pathName'
// 'stat' is also the number of minEntries for MakeResPlotsQTree
//
- if (fDebugLevel > 0) Info("MakeReport","Writing plots and trees to '%s'.", pathName);
+ if (GetDebugLevel() > 0) Info("MakeReport","Writing plots and trees to '%s'.", pathName);
MakeAmpPlots(stat, pathName);
MakeDeltaPlots(pathName);
FitResolutionNew(pathName);
}
-void AliTPCcalibTracks::MakeAmpPlots(Int_t stat, char* pathName){
+void AliTPCcalibTracks::MakeAmpPlots(Int_t stat, const char* pathName){
//
// creates several plots:
// fArrayAmp.ps, fArrayAmpRow.ps and DeltaYZ.ps
allAmpHisOROC->SetTitle("Amp all OROCs");
ps = new TPostScript("fArrayAmp.ps", 112);
- if (fDebugLevel > -1) cout << "creating fArrayAmp.ps..." << endl;
+ if (GetDebugLevel() > -1) cout << "creating fArrayAmp.ps..." << endl;
for (Int_t i = 0; i < fArrayAmp->GetEntriesFast(); i++){
if ( ((TH1F*)fArrayAmp->At(i))->GetEntries() < stat ) continue;
ps->NewPage();
Double_t min = 0;
Double_t max = 0;
ps = new TPostScript("fArrayAmpRow.ps", 112);
- if (fDebugLevel > -1) cout << "creating fArrayAmpRow.ps..." << endl;
+ if (GetDebugLevel() > -1) cout << "creating fArrayAmpRow.ps..." << endl;
for (Int_t i = 0; i < fArrayAmpRow->GetEntriesFast(); i++){
his = (TH1F*)fArrayAmpRow->At(i);
if (his->GetEntries() < stat) continue;
}
-void AliTPCcalibTracks::MakeDeltaPlots(char* pathName){
+void AliTPCcalibTracks::MakeDeltaPlots(const char* pathName){
//
// creates several plots:
// DeltaYZ.ps: DeltaY and DeltaZ histogram with gaus fit
Double_t max = 0;
ps = new TPostScript("DeltaYZ.ps", 112);
- if (fDebugLevel > -1) cout << "creating DeltaYZ.ps..." << endl;
+ if (GetDebugLevel() > -1) cout << "creating DeltaYZ.ps..." << endl;
min = fDeltaY->GetBinCenter(fDeltaY->GetMaximumBin())-20;
max = fDeltaY->GetBinCenter(fDeltaY->GetMaximumBin())+20;
fDeltaY->SetAxisRange(min, max);
}
-void AliTPCcalibTracks::MakeChargeVsDriftLengthPlotsOld(char* pathName){
+void AliTPCcalibTracks::MakeChargeVsDriftLengthPlotsOld(const char* pathName){
//
// creates charge vs. driftlength plots, one TProfile for each ROC
// is not correct like this, should be one TProfile for each sector and padsize
TCanvas* c1 = new TCanvas(); // valgrind 3 ??? 634 bytes in 28 blocks are still reachable
TPostScript *ps;
ps = new TPostScript("chargeVsDriftlengthOld.ps", 112);
- if (fDebugLevel > -1) cout << "creating chargeVsDriftlength.ps..." << endl;
+ if (GetDebugLevel() > -1) cout << "creating chargeVsDriftlength.ps..." << endl;
TProfile *chargeVsDriftlengthAllIROCs = ((TProfile*)fArrayChargeVsDriftlength->At(0)->Clone());
TProfile *chargeVsDriftlengthAllOROCs = ((TProfile*)fArrayChargeVsDriftlength->At(36)->Clone());
chargeVsDriftlengthAllIROCs->SetName("allAmpHisIROC");
}
-void AliTPCcalibTracks::MakeChargeVsDriftLengthPlots(char* pathName){
+void AliTPCcalibTracks::MakeChargeVsDriftLengthPlots(const char* pathName){
//
// creates charge vs. driftlength plots, one TProfile for each ROC
// under development....
c1->Divide(0,3);
TPostScript *ps;
ps = new TPostScript("chargeVsDriftlength.ps", 111);
- if (fDebugLevel > -1) cout << "creating chargeVsDriftlengthNew.ps..." << endl;
+ if (GetDebugLevel() > -1) cout << "creating chargeVsDriftlengthNew.ps..." << endl;
TProfile *chargeVsDriftlengthAllShortPads = ((TProfile*)fcalPadRegionChargeVsDriftlength->GetObject(0,0)->Clone());
TProfile *chargeVsDriftlengthAllMediumPads = ((TProfile*)fcalPadRegionChargeVsDriftlength->GetObject(0,1)->Clone());
-void AliTPCcalibTracks::FitResolutionNew(char* pathName){
+void AliTPCcalibTracks::FitResolutionNew(const char* pathName){
//
// calculates different resulution fits in Y and Z direction
// the histograms are written to 'ResolutionYZ.ps'
TCanvas c;
c.Divide(2,1);
- if (fDebugLevel > -1) cout << "creating ResolutionYZ.ps..." << endl;
+ if (GetDebugLevel() > -1) cout << "creating ResolutionYZ.ps..." << endl;
TPostScript *ps = new TPostScript("ResolutionYZ.ps", 112);
TF2 *fres = new TF2("fres","TMath::Sqrt([0]*[0]+[1]*[1]*x+[2]*[2]*y*y)",0,250,0,1);
fres->SetParameter(0,0.02);
}
-void AliTPCcalibTracks::FitRMSNew(char* pathName){
+void AliTPCcalibTracks::FitRMSNew(const char* pathName){
//
// calculates different resulution-rms fits in Y and Z direction
// the histograms are written to 'RMS_YZ.ps'
TCanvas c; // valgrind 3 42,120 bytes in 405 blocks are still reachable 23,816 bytes in 229 blocks are still reachable
c.Divide(2,1);
- if (fDebugLevel > -1) cout << "creating RMS_YZ.ps..." << endl;
+ if (GetDebugLevel() > -1) cout << "creating RMS_YZ.ps..." << endl;
TPostScript *ps = new TPostScript("RMS_YZ.ps", 112);
TF2 *frms = new TF2("fres","TMath::Sqrt([0]*[0]+[1]*[1]*x+[2]*[2]*y*y)",0,250,0,1);
frms->SetParameter(0,0.02);
}
-void AliTPCcalibTracks::MakeResPlotsQTree(Int_t minEntries, char* pathName){
+void AliTPCcalibTracks::MakeResPlotsQTree(Int_t minEntries, const char* pathName){
//
// Make tree with resolution parameters
// the result is written to 'resol.root' in directory 'pathname'
// RMSe1: error of sigma of gaus fit in RMS-Histogram
//
- if (fDebugLevel > -1) cout << " ***** this is MakeResPlotsQTree *****" << endl;
- if (fDebugLevel > -1) cout << " relax, the calculation will take a while..." << endl;
+ if (GetDebugLevel() > -1) cout << " ***** this is MakeResPlotsQTree *****" << endl;
+ if (GetDebugLevel() > -1) cout << " relax, the calculation will take a while..." << endl;
gSystem->MakeDirectory(pathName);
gSystem->ChangeDirectory(pathName);
}
}
}
- if (fDebugLevel > -1) cout << "Histograms loaded, starting to proces..." << endl;
+ if (GetDebugLevel() > -1) cout << "Histograms loaded, starting to proces..." << endl;
//--------------------------------------------------------------------------------------------
// loop pad type
for (Int_t iq = -1; iq < 10; iq++){
// LOOP Q
- if (fDebugLevel > -1)
+ if (GetDebugLevel() > -1)
cout << "Loop-counter, this is loop " << loopCounter << " of 66, ("
<< (Int_t)((loopCounter)/66.*100) << "% done), "
<< "idim = " << idim << ", ipad = " << ipad << ", iq = " << iq << " \r" << std::flush;
"RMSe0="<<errorRMS<< // error of mean of gaus fit in RMS-Histogram
"RMSe1="<<errorSigma<< // error of sigma of gaus fit in RMS-Histogram
"\n";
- if (fDebugLevel > 5) {
+ if (GetDebugLevel() > 5) {
projectionRes->SetDirectory(fTreeResol.GetFile());
projectionRes->Write(projectionRes->GetName());
projectionRes->SetDirectory(0);
fileInfo.Write("fileInfo");
// resolFile.Close();
// fTreeResol.GetFile()->Close();
- if (fDebugLevel > -1) cout << endl;
- if (fDebugLevel > -1) cout << "MakeResPlotsQTree done, results are in '"<< kFileName.Data() <<"'." << endl;
- gSystem->ChangeDirectory("..");
-}
-
-
-
-
-
-/*
-Int_t AliTPCcalibTracks::fgLoopCounter = 0;
-void AliTPCcalibTracks::MakeResPlotsQTreeThread(Int_t minEntries, char* pathName){
- //
- //
- //
- if (fDebugLevel > -1) cout << " ***** this is MakeResPlotsQTreeThread *****" << endl;
- if (fDebugLevel > -1) cout << " relax, the calculation will take a while..." << endl;
- if (fDebugLevel > -1) cout << " it will be done using 6 TThreads." << endl;
-
- gSystem->MakeDirectory(pathName);
- gSystem->ChangeDirectory(pathName);
- TString kFileName = "resol.root";
-// TTreeSRedirector *fTreeResol = new TTreeSRedirector(kFileName.Data());
- TTreeSRedirector fTreeResol(kFileName.Data());
-
- TH3F *resArray[2][3][11];
- TH3F *rmsArray[2][3][11];
-
- // load histograms from fArrayQDY and fArrayQDZ
- // into resArray and rmsArray
- // that is all we need here
- for (Int_t idim = 0; idim < 2; idim++){
- for (Int_t ipad = 0; ipad < 3; ipad++){
- for (Int_t iq = 0; iq <= 10; iq++){
- rmsArray[idim][ipad][iq]=0;
- resArray[idim][ipad][iq]=0;
- Int_t bin = GetBin(iq,ipad);
- TH3F *hresl = 0;
- if (idim == 0) hresl = (TH3F*)fArrayQDY->At(bin);
- if (idim == 1) hresl = (TH3F*)fArrayQDZ->At(bin);
- if (!hresl) continue;
- resArray[idim][ipad][iq] = (TH3F*) hresl->Clone();
- resArray[idim][ipad][iq]->SetDirectory(0);
- TH3F * hreslRMS = 0;
- if (idim == 0) hreslRMS = (TH3F*)fArrayQRMSY->At(bin);
- if (idim == 1) hreslRMS = (TH3F*)fArrayQRMSZ->At(bin);
- if (!hreslRMS) continue;
- rmsArray[idim][ipad][iq] = (TH3F*) hreslRMS->Clone();
- rmsArray[idim][ipad][iq]->SetDirectory(0);
- }
- }
- }
- if (fDebugLevel > 4) cout << "Histograms loaded, starting to proces..." << endl;
-
- //--------------------------------------------------------------------------------------------
-
- Int_t threadCounter = 0;
- TObjArray *listOfThreads = new TObjArray();
-
- fgLoopCounter = 0;
- for (Int_t idim = 0; idim < 2; idim++){
- // Loop y-z corrdinate
- for (Int_t ipad = 0; ipad < 3; ipad++){
- // loop pad type
-
- // make list of variables for threads
- // make threads
- // list them
- // execute them
-
- TthreadParameterStruct *structOfParameters = new TthreadParameterStruct();
- structOfParameters->logLevel = fDebugLevel;
- structOfParameters->minEntries = minEntries;
- structOfParameters->dim = idim;
- structOfParameters->pad = ipad;
- structOfParameters->resArray = &resArray;
- structOfParameters->rmsArray = &rmsArray;
- structOfParameters->fileName = &kFileName;
- structOfParameters->fTreeResol = &fTreeResol;
- TThread *thread = new TThread(Form("thread%i", threadCounter), (void(*) (void *))&MakeResPlotsQTreeThreadFunction, (void*)structOfParameters);
- listOfThreads->AddAt(thread, threadCounter);
- thread->Run();
-// gSystem->Sleep(500); // so that the threads do not run synchron
-
-// typedef TH3F test;
-// test *testArray;
-// TH3F* (*testArray)[2][3][11];
-// testArray = &resArray;
- int i[2][3][4];
- int (*ptr)[2][3][4];
- ptr = &i;
- int (*ptr2)[2][3][4] = ptr;
- int j = (*ptr2)[1][1][1];
- j = j;
-
- threadCounter++;
- } // ipad-loop
- } // idim-loop
-
- // wait untill all threads are finished
- Bool_t allFinished = kFALSE;
- Int_t numberOfRunningThreads = 0;
- char c[4] = {'-', '\\', '|', '/'};
- Int_t iTime = 0;
- while (!allFinished) {
- allFinished = kTRUE;
- numberOfRunningThreads = 0;
- for (Int_t i = 0; i < listOfThreads->GetEntriesFast(); i++) {
- if (listOfThreads->At(i) != 0x0 && ((TThread*)listOfThreads->At(i))->GetState() == TThread::kRunningState) {
- allFinished = kFALSE;
- numberOfRunningThreads++;
- }
- }
- cout << "Loop-counter, loop " << fgLoopCounter << " of 66 has just started, ("
- << (Int_t)((fgLoopCounter)/66.*100) << "% done), " << "number of running TThreads: "
- << numberOfRunningThreads << " \t" << c[iTime%4] << " \r" << std::flush;
- iTime++;
- gSystem->Sleep(500);
- }
- cout << endl;
-
- // old version:
- // Real time 0:01:31, CP time 44.690
-
- // version without sleep:
- // Real time 0:02:18, CP time 106.280
-
- // version with sleep, listOfThreads-Bug corrected:
- // Real time 0:01:35, CP time 0.800
-
- TObjString fileInfo(Form("Resolution tree, minEntries = %i", minEntries));
- fileInfo.Write("fileInfo");
- if (fDebugLevel > -1) cout << endl;
- if (fDebugLevel > -1) cout << "MakeResPlotsQTree done, results are in '"<< kFileName.Data() <<"'." << endl;
+ if (GetDebugLevel() > -1) cout << endl;
+ if (GetDebugLevel() > -1) cout << "MakeResPlotsQTree done, results are in '"<< kFileName.Data() <<"'." << endl;
gSystem->ChangeDirectory("..");
}
-TMutex* AliTPCcalibTracks::fgWriteMutex = new TMutex();
-TMutex* AliTPCcalibTracks::fgFitResMutex = new TMutex();
-TMutex* AliTPCcalibTracks::fgFitRmsMutex = new TMutex();
-void* AliTPCcalibTracks::MakeResPlotsQTreeThreadFunction(void* arg){
- //
- //
- //
-
- TthreadParameterStruct *structOfParameters = (TthreadParameterStruct*)arg;
- Int_t fDebugLevel = structOfParameters->logLevel;
- Int_t minEntries = structOfParameters->minEntries;
- Int_t idim = structOfParameters->dim;
- Int_t ipad = structOfParameters->pad;
- TThread::Lock();
- TH3F* (*resArray)[2][3][11] = structOfParameters->resArray;
- TH3F* (*rmsArray)[2][3][11] = structOfParameters->rmsArray;
- TThread::UnLock();
- TString *kFileName = structOfParameters->fileName;
- TTreeSRedirector *fTreeResol = structOfParameters->fTreeResol;
-
- if (fDebugLevel > 4) TThread::Printf("Thread started, dim = %i, pad = %i...", idim, ipad);
-
- TThread::Lock();
- char name[200];
- sprintf(name, "dim%ipad%i", idim, ipad);
- TH1D *projectionRes = new TH1D(Form("projectionRes%s", name), "projectionRes", 50, -1, 1);
- TH1D *projectionRms = new TH1D(Form("projectionRms%s", name), "projectionRms", 50, -1, 1);
- char fitFuncName[200];
- sprintf(name, "fitFunctionDim%iPad%i", idim, ipad);
- TF1 *fitFunction = new TF1(fitFuncName, "gaus");
- TThread::UnLock();
-
- Double_t zMean, angleMean, zCenter, angleCenter;
- Double_t zSigma, angleSigma;
-
- for (Int_t iq = -1; iq < 10; iq++){
- // LOOP Q
- fgLoopCounter++;
- TH3F *hres = 0;
- TH3F *hrms = 0;
- Double_t qMean = 0;
- Float_t entriesQ = 0;
-
- if (fDebugLevel > 4) TThread::Printf(" start of iq-loop, dim = %i, pad = %i, iq = %i...", idim, ipad, iq);
- // calculate qMean
- if (iq == -1){
- // integrated spectra
- for (Int_t iql = 0; iql < 10; iql++){
- Int_t bin = GetBin(iql,ipad);
- TH3F *hresl = (*resArray)[idim][ipad][iql];
- TH3F *hrmsl = (*rmsArray)[idim][ipad][iql];
- if (!hresl) continue;
- if (!hrmsl) continue;
- entriesQ += hresl->GetEntries();
- qMean += hresl->GetEntries() * GetQ(bin);
- if (!hres) {
- TThread::Lock();
- hres = (TH3F*)hresl->Clone();
- hrms = (TH3F*)hrmsl->Clone();
- TThread::UnLock();
- }
- else{
- hres->Add(hresl);
- hrms->Add(hrmsl);
- }
- }
- qMean /= entriesQ;
- qMean *= -1.; // integral mean charge
- }
- else {
- // loop over neighboured Q-bins
- // accumulate entries from neighboured Q-bins
- for (Int_t iql = iq - 1; iql <= iq + 1; iql++){
- if (iql < 0) continue;
- Int_t bin = GetBin(iql,ipad);
- TH3F * hresl = (*resArray)[idim][ipad][iql];
- TH3F * hrmsl = (*rmsArray)[idim][ipad][iql];
- if (!hresl) continue;
- if (!hrmsl) continue;
- entriesQ += hresl->GetEntries();
- qMean += hresl->GetEntries() * GetQ(bin);
- if (!hres) {
- TThread::Lock();
- hres = (TH3F*) hresl->Clone();
- hrms = (TH3F*) hrmsl->Clone();
- TThread::UnLock();
- }
- else{
- hres->Add(hresl);
- hrms->Add(hrmsl);
- }
- }
- qMean/=entriesQ;
- }
- TAxis *xAxisDriftLength = hres->GetXaxis(); // driftlength / z - axis
- TAxis *yAxisAngle = hres->GetYaxis(); // angle axis
- TAxis *zAxisDelta = hres->GetZaxis(); // delta axis
- TAxis *zAxisRms = hrms->GetZaxis(); // rms axis
-
- // loop over all angle bins
- for (Int_t ibinyAngle = 1; ibinyAngle <= yAxisAngle->GetNbins(); ibinyAngle++) {
- angleCenter = yAxisAngle->GetBinCenter(ibinyAngle);
- // loop over all driftlength bins
- for (Int_t ibinxDL = 1; ibinxDL <= xAxisDriftLength->GetNbins(); ibinxDL++) {
- zCenter = xAxisDriftLength->GetBinCenter(ibinxDL);
- zSigma = xAxisDriftLength->GetBinWidth(ibinxDL);
- angleSigma = yAxisAngle->GetBinWidth(ibinyAngle);
- zMean = zCenter; // changens, when more statistic is accumulated
- angleMean = angleCenter; // changens, when more statistic is accumulated
-
- // create 2 1D-Histograms, projectionRes and projectionRms
- // these histograms are delta histograms for given direction, padSize, chargeBin,
- // angleBin and driftLengthBin
- // later on they will be fitted with a gausian, its sigma is the resoltuion...
- sprintf(name,"%s, zCenter: %f, angleCenter: %f", hres->GetName(), zCenter, angleCenter);
- projectionRes->SetNameTitle(name, name);
- sprintf(name,"%s, zCenter: %f, angleCenter: %f", hrms->GetName(),zCenter,angleCenter);
- projectionRms->SetNameTitle(name, name);
-
- projectionRes->Reset();
- projectionRes->SetBins(zAxisDelta->GetNbins(), zAxisDelta->GetXmin(), zAxisDelta->GetXmax());
- projectionRms->Reset();
- projectionRms->SetBins(zAxisRms->GetNbins(), zAxisRms->GetXmin(), zAxisRms->GetXmax());
- projectionRes->SetDirectory(0);
- projectionRms->SetDirectory(0);
-
- Double_t entries = 0;
- Int_t nbins = 0; // counts, how many bins were accumulated
- zMean = 0;
- angleMean = 0;
- entries = 0;
-
- // fill projectionRes and projectionRms for given dim, ipad and iq,
- // as well as for given angleBin and driftlengthBin
- // if this gives not enough statistic, include neighbourhood
- // (angle and driftlength) successifely
- for (Int_t dbin = 0; dbin <= 8; dbin++){ // delta-bins around centered angleBin and driftlengthBin
- for (Int_t dbiny2 = -1; dbiny2 <= 1; dbiny2++) { // delta-bins in angle direction
- for (Int_t dbinx2 = -3; dbinx2 <= 3; dbinx2++){ // delta-bins in driftlength direction
- if (TMath::Abs(dbinx2) + TMath::Abs(dbiny2) != dbin) continue; // add each bin only one time !
- Int_t binx2 = ibinxDL + dbinx2; // position variable in x (driftlength) direction
- Int_t biny2 = ibinyAngle + dbiny2; // position variable in y (angle) direction
- if (binx2 < 1 || biny2 < 1) continue; // don't go out of the histogram!
- if (binx2 >= xAxisDriftLength->GetNbins()) continue; // don't go out of the histogram!
- if (biny2 >= yAxisAngle->GetNbins()) continue; // don't go out of the histogram!
- nbins++; // count the number of accumulated bins
- // Fill resolution histo
- for (Int_t ibin3 = 1; ibin3 < zAxisDelta->GetNbins(); ibin3++) {
- // Int_t content = (Int_t)hres->GetBinContent(binx2, biny2, ibin3); // unused variable
- projectionRes->Fill(zAxisDelta->GetBinCenter(ibin3), hres->GetBinContent(binx2, biny2, ibin3));
- entries += hres->GetBinContent(binx2, biny2, ibin3);
- zMean += hres->GetBinContent(binx2, biny2, ibin3) * xAxisDriftLength->GetBinCenter(binx2);
- angleMean += hres->GetBinContent(binx2, biny2, ibin3) * yAxisAngle->GetBinCenter(biny2);
- } // ibin3 loop
- // fill RMS histo
- for (Int_t ibin3 = 1; ibin3 < zAxisRms->GetNbins(); ibin3++) {
- projectionRms->Fill(zAxisRms->GetBinCenter(ibin3), hrms->GetBinContent(binx2, biny2, ibin3));
- }
- } //dbinx2 loop
- if (entries > minEntries) break; // enough statistic accumulated
- } // dbiny2 loop
- if (entries > minEntries) break; // enough statistic accumulated
- } // dbin loop
- if ( entries< minEntries) continue; // when it was absolutly impossible to get enough statistic, don't write this point into the resolution tree
- zMean /= entries;
- angleMean /= entries;
-
- if (entries > minEntries) {
- // when enough statistic is accumulated
- // fit Delta histograms with a gausian
- // of the gausian is the resolution (resol), its fit error is sigma
- // store also mean and RMS of the histogram
- Float_t xmin = projectionRes->GetMean() - 2. * projectionRes->GetRMS() - 0.2;
- Float_t xmax = projectionRes->GetMean() + 2. * projectionRes->GetRMS() + 0.2;
- fitFunction->SetMaximum(xmax);
- fitFunction->SetMinimum(xmin);
- TThread::Lock();
- // fgFitResMutex->Lock();
- projectionRes->Fit(fitFuncName, "qN0", "", xmin, xmax);
- // fgFitResMutex->UnLock();
- TThread::UnLock();
- Float_t resol = fitFunction->GetParameter(2);
- Float_t sigma = fitFunction->GetParError(2);
- Float_t meanR = projectionRes->GetMean();
- Float_t sigmaR = projectionRes->GetRMS();
- // fit also RMS histograms with a gausian
- // store mean and sigma of the gausian in rmsMean and rmsSigma
- // store also the fit errors in errorRMS and errorSigma
- xmin = projectionRms->GetMean() - 2. * projectionRes->GetRMS() - 0.2;
- xmax = projectionRms->GetMean() + 2. * projectionRes->GetRMS() + 0.2;
- fitFunction->SetMaximum(xmax);
- fitFunction->SetMinimum(xmin);
- TThread::Lock();
- projectionRms->Fit(fitFuncName, "qN0", "", xmin, xmax);
- TThread::UnLock();
- Float_t rmsMean = fitFunction->GetParameter(1);
- Float_t rmsSigma = fitFunction->GetParameter(2);
- Float_t errorRMS = fitFunction->GetParError(1);
- Float_t errorSigma = fitFunction->GetParError(2);
- Float_t length = 0.75;
- if (ipad == 1) length = 1;
- if (ipad == 2) length = 1.5;
-
- fgWriteMutex->Lock();
- (*fTreeResol)<<"Resol"<<
- "Entries="<<entries<< // number of entries for this resolution point
- "nbins="<<nbins<< // number of bins that were accumulated
- "Dim="<<idim<< // direction, Dim==0: y-direction, Dim==1: z-direction
- "Pad="<<ipad<< // padSize; short, medium and long
- "Length="<<length<< // pad length, 0.75, 1, 1.5
- "QMean="<<qMean<< // mean charge of current charge bin and its neighbours, Qmean<0: integrated spectra
- "Zc="<<zCenter<< // center of middle bin in drift direction
- "Zm="<<zMean<< // mean dirftlength for accumulated Delta-Histograms
- "Zs="<<zSigma<< // width of driftlength bin
- "AngleC="<<angleCenter<< // center of middle bin in Angle-Direction
- "AngleM="<<angleMean<< // mean angle for accumulated Delta-Histograms
- "AngleS="<<angleSigma<< // width of Angle-bin
- "Resol="<<resol<< // sigma for gaus fit through Delta-Histograms
- "Sigma="<<sigma<< // error of sigma for gaus fit through Delta Histograms
- "MeanR="<<meanR<< // mean of the Delta-Histogram
- "SigmaR="<<sigmaR<< // rms of the Delta-Histogram
- "RMSm="<<rmsMean<< // mean of the gaus fit through RMS-Histogram
- "RMSs="<<rmsSigma<< // sigma of the gaus fit through RMS-Histogram
- "RMSe0="<<errorRMS<< // error of mean of gaus fit in RMS-Histogram
- "RMSe1="<<errorSigma<< // error of sigma of gaus fit in RMS-Histogram
- "\n";
- if (fDebugLevel > 5) {
- projectionRes->SetDirectory(fTreeResol->GetFile());
- projectionRes->Write(projectionRes->GetName());
- projectionRes->SetDirectory(0);
- projectionRms->SetDirectory(fTreeResol->GetFile());
- projectionRms->Write(projectionRms->GetName());
- projectionRes->SetDirectory(0);
- }
- fgWriteMutex->UnLock();
- } // if (projectionRes->GetSum() > minEntries)
- } // for (Int_t ibinxDL = 1; ibinxDL <= xAxisDriftLength->GetNbins(); ibinxDL++)
- } // for (Int_t ibinyAngle = 1; ibinyAngle <= yAxisAngle->GetNbins(); ibinyAngle++)
-
- } // iq-loop
- delete projectionRes;
- delete projectionRms;
- if (fDebugLevel > 4) TThread::Printf("Ende, dim = %i, pad = %i", idim, ipad);
- return 0;
-}
-*/
Long64_t AliTPCcalibTracks::Merge(TCollection *collectionList) {
// Be carefull: histograms are linked to a file, switch this off by TH1::AddDirectory(kFALSE) !!!
//
- if (fDebugLevel > 0) cout << " ***** this is AliTPCcalibTracks::Merge(TCollection *collectionList) *****"<< endl;
+ if (GetDebugLevel() > 0) cout << " ***** this is AliTPCcalibTracks::Merge(TCollection *collectionList) *****"<< endl;
if (!collectionList) return 0;
if (collectionList->IsEmpty()) return -1;
- if (fDebugLevel > 1) cout << "the collectionList contains " << collectionList->GetEntries() << " entries." << endl; // REMOVE THIS LINE!!!!!!!!!!!!!!!!!1
- if (fDebugLevel > 5) cout << " the list in the merge-function looks as follows: " << endl;
+ if (GetDebugLevel() > 1) cout << "the collectionList contains " << collectionList->GetEntries() << " entries." << endl; // REMOVE THIS LINE!!!!!!!!!!!!!!!!!1
+ if (GetDebugLevel() > 5) cout << " the list in the merge-function looks as follows: " << endl;
collectionList->Print();
// create a list for each data member
TIterator *listIterator = collectionList->MakeIterator();
AliTPCcalibTracks *calibTracks = 0;
- if (fDebugLevel > 1) cout << "start to iterate, filling lists" << endl;
+ if (GetDebugLevel() > 1) cout << "start to iterate, filling lists" << endl;
Int_t counter = 0;
- while ( (calibTracks = (AliTPCcalibTracks*)listIterator->Next()) ){
+ while ( (calibTracks = dynamic_cast<AliTPCcalibTracks*> (listIterator->Next())) ){
// loop over all entries in the collectionList and get dataMembers into lists
if (!calibTracks) continue;
clusterCutHistoList->Add(calibTracks->GetfClusterCutHisto());
hclusterPerPadrowList->Add(calibTracks->GetfHclusterPerPadrow());
hclusterPerPadrowRawList->Add(calibTracks->GetfHclusterPerPadrowRaw());
- fCalPadClusterPerPad->Add(calibTracks->GetfCalPadClusterPerPad());
- fCalPadClusterPerPadRaw->Add(calibTracks->GetfCalPadClusterPerPadRaw());
+ //
+ if (fCalPadClusterPerPad && calibTracks->GetfCalPadClusterPerPad())
+ fCalPadClusterPerPad->Add(calibTracks->GetfCalPadClusterPerPad());
+ // fCalPadClusterPerPadRaw->Add(calibTracks->GetfCalPadClusterPerPadRaw());
counter++;
- if (fDebugLevel > 5) cout << "filling lists, object " << counter << " added." << endl;
+ if (GetDebugLevel() > 5) cout << "filling lists, object " << counter << " added." << endl;
}
// merge data members
- if (fDebugLevel > 0) cout << "histogram's merge-functins are called... " << endl;
+ if (GetDebugLevel() > 0) cout << "histogram's merge-functins are called... " << endl;
fDeltaY->Merge(deltaYList);
fDeltaZ->Merge(deltaZList);
fHclus->Merge(hclusList);
TList* histList = 0;
TIterator *objListIterator = 0;
- if (fDebugLevel > 0) cout << "merging fArrayAmpRows..." << endl;
+ if (GetDebugLevel() > 0) cout << "merging fArrayAmpRows..." << endl;
// merge fArrayAmpRows
for (Int_t i = 0; i < fArrayAmpRow->GetEntriesFast(); i++ ) { // loop over data member, i<72
objListIterator = arrayAmpRowList->MakeIterator();
delete objListIterator;
}
- if (fDebugLevel > 0) cout << "merging fArrayAmps..." << endl;
+ if (GetDebugLevel() > 0) cout << "merging fArrayAmps..." << endl;
// merge fArrayAmps
for (Int_t i = 0; i < fArrayAmp->GetEntriesFast(); i++ ) { // loop over data member, i<72
- TIterator *objListIterator = arrayAmpList->MakeIterator();
+ TIterator *cobjListIterator = arrayAmpList->MakeIterator();
histList = new TList;
- while (( objarray = (TObjArray*)objListIterator->Next() )) {
+ while (( objarray = (TObjArray*)cobjListIterator->Next() )) {
// loop over arrayAmpList, get TObjArray, get object at position i, cast it into TH1F
if (!objarray) continue;
hist = (TH1F*)(objarray->At(i));
}
((TH1F*)(fArrayAmp->At(i)))->Merge(histList);
delete histList;
- delete objListIterator;
+ delete cobjListIterator;
}
- if (fDebugLevel > 0) cout << "merging fArrayQDY..." << endl;
+ if (GetDebugLevel() > 0) cout << "merging fArrayQDY..." << endl;
// merge fArrayQDY
for (Int_t i = 0; i < fArrayQDY->GetEntriesFast(); i++) { // loop over data member, i < 300
objListIterator = arrayQDYList->MakeIterator();
delete objListIterator;
}
- if (fDebugLevel > 0) cout << "merging fArrayQDZ..." << endl;
+ if (GetDebugLevel() > 0) cout << "merging fArrayQDZ..." << endl;
// merge fArrayQDZ
for (Int_t i = 0; i < fArrayQDZ->GetEntriesFast(); i++) { // loop over data member, i < 300
objListIterator = arrayQDZList->MakeIterator();
delete objListIterator;
}
- if (fDebugLevel > 0) cout << "merging fArrayQRMSY..." << endl;
+ if (GetDebugLevel() > 0) cout << "merging fArrayQRMSY..." << endl;
// merge fArrayQRMSY
for (Int_t i = 0; i < fArrayQRMSY->GetEntriesFast(); i++) { // loop over data member, i < 300
objListIterator = arrayQRMSYList->MakeIterator();
delete objListIterator;
}
- if (fDebugLevel > 0) cout << "merging fArrayQRMSZ..." << endl;
+ if (GetDebugLevel() > 0) cout << "merging fArrayQRMSZ..." << endl;
// merge fArrayQRMSZ
for (Int_t i = 0; i < fArrayQRMSZ->GetEntriesFast(); i++) { // loop over data member, i < 300
objListIterator = arrayQRMSZList->MakeIterator();
delete objListIterator;
}
- if (fDebugLevel > 0) cout << "merging fArrayChargeVsDriftlength..." << endl;
+ if (GetDebugLevel() > 0) cout << "merging fArrayChargeVsDriftlength..." << endl;
// merge fArrayChargeVsDriftlength
for (Int_t i = 0; i < fArrayChargeVsDriftlength->GetEntriesFast(); i++) { // loop over data member, i < 300
objListIterator = arrayChargeVsDriftlengthList->MakeIterator();
delete objListIterator;
}
- if (fDebugLevel > 0) cout << "merging fcalPadRegionChargeVsDriftlength..." << endl;
+ if (GetDebugLevel() > 0) cout << "merging fcalPadRegionChargeVsDriftlength..." << endl;
// merge fcalPadRegionChargeVsDriftlength
AliTPCCalPadRegion *cpr = 0x0;
- if (fDebugLevel > 0) cout << "starting to merge the rest: fResolY, fResolZ , fRMSY, fRMSZ..." << endl;
+ if (GetDebugLevel() > 0) cout << "starting to merge the rest: fResolY, fResolZ , fRMSY, fRMSZ..." << endl;
// merge fResolY
for (Int_t i = 0; i < fResolY->GetEntriesFast(); i++) { // loop over data member, i < 3
objListIterator = resolYList->MakeIterator();
delete resolZList;
delete rMSYList;
delete rMSZList;
-// delete nRowsList;
-// delete nSectList;
-// delete fileNoList;
delete listIterator;
- if (fDebugLevel > 0) cout << "merging done!" << endl;
+ if (GetDebugLevel() > 0) cout << "merging done!" << endl;
return 1;
}
}
+void AliTPCcalibTracks::MakeQPosNormAll(TTree * chainres, AliTPCClusterParam * param, Int_t maxPoints){
+ //
+ // Make position corrections
+ // for the moment Only using debug streamer
+ // chainres - debug tree
+ // param - parameters to be updated
+ // maxPoints - maximal number of points using for fit
+ // verbose - print info flag
+ //
+ // Current implementation - only using debug streamers
+ //
+
+ /*
+ //Defaults
+ Int_t maxPoints=100000;
+ */
+ /*
+ Usage:
+ //0. Load libraries
+ gSystem->Load("libANALYSIS");
+ gSystem->Load("libSTAT");
+ gSystem->Load("libTPCcalib");
+
+
+ //1. Load Parameters to be modified:
+ //e.g:
+
+ AliCDBManager::Instance()->SetDefaultStorage("local://$ALICE_ROOT");
+ AliCDBManager::Instance()->SetRun(0)
+ AliTPCClusterParam * param = AliTPCcalibDB::Instance()->GetClusterParam();
+
+ //2. Load chain from debug streamers
+ //
+ //e.g
+ gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
+ gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
+ AliXRDPROOFtoolkit tool;
+ TChain * chainres = tool.MakeChain("tracks.txt","ResolCl",0,10200);
+ chainres->Lookup();
+ //3. Do fits and store results
+ //
+ AliTPCcalibTracks::MakeQPosNormAll(chainres,param,200000,0) ;
+ TFile f("paramout.root","recreate");
+ param->Write("clusterParam");
+ f.Close();
+ //4. Verification
+ TFile f2("paramout.root");
+ AliTPCClusterParam *param2 = (AliTPCClusterParam*)f2.Get("clusterParam");
+ param2->SetInstance(param2);
+ chainres->Draw("fitZ0:AliTPCClusterParam::SPosCorrection(1,0,Cl.fPad,Cl.fTimeBin,Cl.fZ,Cl.fSigmaY2,Cl.fSigmaZ2,Cl.fMax)","Cl.fDetector<36","",10000); // should be line
+
+ */
+
+
+ TStatToolkit toolkit;
+ Double_t chi2;
+ TVectorD fitParamY0;
+ TVectorD fitParamY1;
+ TVectorD fitParamZ0;
+ TVectorD fitParamZ1;
+ TMatrixD covMatrix;
+ Int_t npoints;
+
+ chainres->SetAlias("dp","(-1+(Cl.fZ>0)*2)*((Cl.fPad-int(Cl.fPad))-0.5)");
+ chainres->SetAlias("dt","(-1+(Cl.fZ>0)*2)*((Cl.fTimeBin-0.66-int(Cl.fTimeBin-0.66))-0.5)");
+ chainres->SetAlias("sp","(sin(dp*pi)-dp*pi)");
+ chainres->SetAlias("st","(sin(dt)-dt)");
+ //
+ chainres->SetAlias("di","sqrt(1.-abs(Cl.fZ/250.))");
+ //
+ //
+ //
+ TCut cutA("1");
+ TString fstringY="";
+ //
+ fstringY+="(dp)++"; //1
+ fstringY+="(dp)*di++"; //2
+ fstringY+="(sp)++"; //3
+ fstringY+="(sp)*di++"; //4
+ TString fstringZ="";
+ fstringZ+="(dt)++"; //1
+ fstringZ+="(dt)*di++"; //2
+ fstringZ+="(st)++"; //3
+ fstringZ+="(st)*di++"; //4
+ //
+ // Z corrections
+ //
+ TString *strZ0 = toolkit.FitPlane(chainres,"(Cl.fZ-PZ0.fElements[0]):CSigmaZ0",fstringZ.Data(), "Cl.fDetector<36"+cutA, chi2,npoints,fitParamZ0,covMatrix,-1,0,maxPoints);
+ printf("Z0 - chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
+ param->fPosZcor[0] = (TVectorD*) fitParamZ0.Clone();
+ //
+ TString *strZ1 = toolkit.FitPlane(chainres,"(Cl.fZ-PZ0.fElements[0]):CSigmaZ0",fstringZ.Data(), "Cl.fDetector>36"+cutA, chi2,npoints,fitParamZ1,covMatrix,-1,0,maxPoints);
+ //
+ printf("Z1 - chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
+ param->fPosZcor[1] = (TVectorD*) fitParamZ1.Clone();
+ param->fPosZcor[2] = (TVectorD*) fitParamZ1.Clone();
+ //
+ // Y corrections
+ //
+ TString *strY0 = toolkit.FitPlane(chainres,"(Cl.fY-PY0.fElements[0]):CSigmaY0",fstringY.Data(), "Cl.fDetector<36"+cutA, chi2,npoints,fitParamY0,covMatrix,-1,0,maxPoints);
+ printf("Y0 - chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
+ param->fPosYcor[0] = (TVectorD*) fitParamY0.Clone();
+
+
+ TString *strY1 = toolkit.FitPlane(chainres,"(Cl.fY-PY0.fElements[0]):CSigmaY0",fstringY.Data(), "Cl.fDetector>36"+cutA, chi2,npoints,fitParamY1,covMatrix,-1,0,maxPoints);
+ //
+ printf("Y1 - chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
+ param->fPosYcor[1] = (TVectorD*) fitParamY1.Clone();
+ param->fPosYcor[2] = (TVectorD*) fitParamY1.Clone();
+ //
+ //
+ //
+ chainres->SetAlias("fitZ0",strZ0->Data());
+ chainres->SetAlias("fitZ1",strZ1->Data());
+ chainres->SetAlias("fitY0",strY0->Data());
+ chainres->SetAlias("fitY1",strY1->Data());
+ // chainres->Draw("Cl.fZ-PZ0.fElements[0]","CSigmaY0<0.7&&CSigmaZ0<0.7"+cutA,"",10000);
+}
+
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff