X-Git-Url: http://git.uio.no/git/?p=u%2Fmrichter%2FAliRoot.git;a=blobdiff_plain;f=TPC%2FAliTPCcalibTracksGain.cxx;h=4d4de07bd752c650d0096e146f2a036070ca8b76;hp=5ffef1af7dcbca624576f7ee2323a6f3577b4d4f;hb=6fc3dfc08ab2e4943b1b70379ac1dc22e202c29b;hpb=c1418a4c2784245bd8d9b17d8c3d297d51328b91 diff --git a/TPC/AliTPCcalibTracksGain.cxx b/TPC/AliTPCcalibTracksGain.cxx index 5ffef1af7dc..4d4de07bd75 100644 --- a/TPC/AliTPCcalibTracksGain.cxx +++ b/TPC/AliTPCcalibTracksGain.cxx @@ -96,8 +96,79 @@ // the fitters after loading this object from file. // (This will be gone for a new ROOT version > v5-17-05) // +// +// In order to debug some numerical algorithm all data data which are used for +// fitters can be stored in the debug streamers. In case of fitting roblems the +// errors and tendencies can be checked. +// +// Debug Streamers: +// +// +// +// +// //////////////////////////////////////////////////////////////////////////// +/* + .x ~/UliStyle.C + gSystem->Load("libANALYSIS"); + gSystem->Load("libSTAT"); + gSystem->Load("libTPCcalib"); + + TFile fcalib("CalibObjects.root"); + TObjArray * array = (TObjArray*)fcalib.Get("TPCCalib"); + AliTPCcalibTracksGain * gain = ( AliTPCcalibTracksGain *)array->FindObject("calibTracksGain"); + + // + // Angular and drift correction + // + AliTPCClusterParam *param = new AliTPCClusterParam;param->SetInstance(param); + gain->UpdateClusterParam(param); + TF2 fdrty("fdrty","AliTPCClusterParam::SQnorm(0,0,x,y,0)",0,1,0,1) + + // + // Make visual Tree - compare with Kr calibration + // + AliTPCCalPad * gain010 = gain->CreateFitCalPad(0,kTRUE,0); gain010->SetName("CGain010"); + AliTPCCalPad * gain110 = gain->CreateFitCalPad(1,kTRUE,0); gain110->SetName("CGain110"); + AliTPCCalPad * gain210 = gain->CreateFitCalPad(2,kTRUE,0); gain210->SetName("CGain210"); + TFile fkr("/u/miranov/GainMap.root"); + AliTPCCalPad *gainKr = fkr.Get("GainMap"); fkr->SetName("KrGain"); + // + AliTPCPreprocessorOnline * preprocesor = new AliTPCPreprocessorOnline; + preprocesor->AddComponent(gain010); + preprocesor->AddComponent(gain110); + preprocesor->AddComponent(gain210); + preprocesor->AddComponent(gainKr); + preprocesor->DumpToFile("cosmicGain.root"); + // + // + // + // Simple session using the debug streamers + // + + gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros"); + gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+") + AliXRDPROOFtoolkit tool; + + TChain * chain0 = tool.MakeChain("gain.txt","dEdx",0,1000000); + TChain * chain1 = tool.MakeChain("gain.txt","Track",0,1000000); + TChain * chain2 = tool.MakeChain("gain.txt","TrackG",0,1000000); + chain0->Lookup(); + chain1->Lookup(); + chain2->Lookup(); + + chain2->SetAlias("k1","1/0.855"); + chain2->SetAlias("k0","1/0.9928"); + chain2->SetAlias("k2","1/1.152"); + +*/ + + + +#include "AliTPCcalibTracksGain.h" + + #include #include #include "TSystem.h" @@ -108,7 +179,6 @@ #include "AliTPCClusterParam.h" #include "AliTrackPointArray.h" #include "TCint.h" -#include "AliTPCcalibTracksGain.h" #include #include #include @@ -119,6 +189,7 @@ #include #include #include +#include // // AliRoot includes @@ -131,6 +202,7 @@ #include "AliTPCCalROC.h" #include "AliTPCCalPad.h" #include "AliTPCClusterParam.h" +#include "AliTPCcalibDB.h" // #include "AliTracker.h" #include "AliESD.h" @@ -141,8 +213,11 @@ #include "AliTPCclusterMI.h" #include "AliTPCcalibTracksCuts.h" #include "AliTPCFitPad.h" +#include "TStatToolkit.h" +#include "TString.h" +#include "TCut.h" -// REMOVE ALL OF THIS +// #include #include "AliESDEvent.h" @@ -171,35 +246,32 @@ AliTPCcalibTracksGain::AliTPCcalibTracksGain() : AliTPCcalibBase(), fCuts(0), // cuts that are used for sieving the tracks used for calibration // - // Simple Array of histograms - // - fArrayQM(0), // Qmax normalized - fArrayQT(0), // Qtot normalized - fProfileArrayQM(0), // Qmax normalized versus local X - fProfileArrayQT(0), // Qtot normalized versus local X - fProfileArrayQM2D(0), // Qmax normalized versus local X and phi - fProfileArrayQT2D(0), // Qtot normalized versus local X and phi - // // Fitters // fSimpleFitter(0), // simple fitter for short pads fSqrtFitter(0), // sqrt fitter for medium pads fLogFitter(0), // log fitter for long pads + fFitter0M(0), // fitting of the atenuation, angular correction, and mean chamber gain fFitter1M(0), // fitting of the atenuation, angular correction, and mean chamber gain fFitter2M(0), // fitting of the atenuation, angular correction, and mean chamber gain fFitter0T(0), // fitting of the atenuation, angular correction, and mean chamber gain fFitter1T(0), // fitting of the atenuation, angular correction, and mean chamber gain fFitter2T(0), // fitting of the atenuation, angular correction, and mean chamber gain + // + fDFitter0M(0), // fitting of the atenuation, angular correction + fDFitter1M(0), // fitting of the atenuation, angular correction + fDFitter2M(0), // fitting of the atenuation, angular correction + fDFitter0T(0), // fitting of the atenuation, angular correction + fDFitter1T(0), // fitting of the atenuation, angular correction + fDFitter2T(0), // fitting of the atenuation, angular correction + // fSingleSectorFitter(0), // just for debugging // // Counters // fTotalTracks(0), // just for debugging - fAcceptedTracks(0), // just for debugging - fDebugCalPadRaw(0), // just for debugging - fDebugCalPadCorr(0), // just for debugging - fPrevIter(0) // the calibration object in its previous iteration (will not be owned by the new object, don't forget to delete it!) + fAcceptedTracks(0) // just for debugging { // @@ -210,15 +282,6 @@ AliTPCcalibTracksGain::AliTPCcalibTracksGain() : AliTPCcalibTracksGain::AliTPCcalibTracksGain(const AliTPCcalibTracksGain& obj) : AliTPCcalibBase(obj), fCuts(obj.fCuts), // cuts that are used for sieving the tracks used for calibration - fArrayQM(0), // Qmax normalized - fArrayQT(0), // Qtot normalized - // - // Simple Histograms - // - fProfileArrayQM(obj.fProfileArrayQM), // Qmax normalized versus local X - fProfileArrayQT(obj.fProfileArrayQT), // Qtot normalized versus local X - fProfileArrayQM2D(obj.fProfileArrayQM2D), // Qmax normalized versus local X and phi - fProfileArrayQT2D(obj.fProfileArrayQT2D), // Qtot normalized versus local X and phi // // Fitters // @@ -231,15 +294,19 @@ AliTPCcalibTracksGain::AliTPCcalibTracksGain(const AliTPCcalibTracksGain& obj) : fFitter0T(obj.fFitter0T), fFitter1T(obj.fFitter1T), fFitter2T(obj.fFitter2T), + // + fDFitter0M(obj.fDFitter0M), + fDFitter1M(obj.fDFitter1M), + fDFitter2M(obj.fDFitter2M), + fDFitter0T(obj.fDFitter0T), + fDFitter1T(obj.fDFitter1T), + fDFitter2T(obj.fDFitter2T), fSingleSectorFitter(obj.fSingleSectorFitter), // just for debugging // // Counters // fTotalTracks(obj.fTotalTracks), // just for debugging - fAcceptedTracks(obj.fAcceptedTracks), // just for debugging - fDebugCalPadRaw(obj.fDebugCalPadRaw), // just for debugging - fDebugCalPadCorr(obj.fDebugCalPadCorr), // just for debugging - fPrevIter(obj.fPrevIter) // the calibration object in its previous iteration (will not be owned by the new object, don't forget to delete it!) + fAcceptedTracks(obj.fAcceptedTracks) // just for debugging { // @@ -254,30 +321,18 @@ AliTPCcalibTracksGain& AliTPCcalibTracksGain::operator=(const AliTPCcalibTracksG if (this != &rhs) { TNamed::operator=(rhs); - fDebugCalPadRaw = new AliTPCCalPad(*(rhs.fDebugCalPadRaw)); - fDebugCalPadCorr = new AliTPCCalPad(*(rhs.fDebugCalPadCorr)); fSimpleFitter = new AliTPCFitPad(*(rhs.fSimpleFitter)); fSqrtFitter = new AliTPCFitPad(*(rhs.fSqrtFitter)); fLogFitter = new AliTPCFitPad(*(rhs.fLogFitter)); fSingleSectorFitter = new AliTPCFitPad(*(rhs.fSingleSectorFitter)); - fPrevIter = new AliTPCcalibTracksGain(*(rhs.fPrevIter)); fCuts = new AliTPCcalibTracksCuts(*(rhs.fCuts)); } return *this; } -AliTPCcalibTracksGain::AliTPCcalibTracksGain(const char* name, const char* title, AliTPCcalibTracksCuts* cuts, TNamed* /*debugStreamPrefix*/, AliTPCcalibTracksGain* prevIter) : +AliTPCcalibTracksGain::AliTPCcalibTracksGain(const char* name, const char* title, AliTPCcalibTracksCuts* cuts) : AliTPCcalibBase(), fCuts(0), // cuts that are used for sieving the tracks used for calibration - fArrayQM(0), // Qmax normalized - fArrayQT(0), // Qtot normalized - // - // Simple Histograms - // - fProfileArrayQM(0), // Qmax normalized versus local X - fProfileArrayQT(0), // Qtot normalized versus local X - fProfileArrayQM2D(0), // Qmax normalized versus local X and phi - fProfileArrayQT2D(0), // Qtot normalized versus local X and phi // // Fitters // @@ -290,24 +345,26 @@ AliTPCcalibTracksGain::AliTPCcalibTracksGain(const char* name, const char* title fFitter0T(0), // fitting of the atenuation, angular correction, and mean chamber gain fFitter1T(0), // fitting of the atenuation, angular correction, and mean chamber gain fFitter2T(0), // fitting of the atenuation, angular correction, and mean chamber gain + // + fDFitter0M(0), // fitting of the atenuation, angular correction + fDFitter1M(0), // fitting of the atenuation, angular correction + fDFitter2M(0), // fitting of the atenuation, angular correction + fDFitter0T(0), // fitting of the atenuation, angular correction + fDFitter1T(0), // fitting of the atenuation, angular correction + fDFitter2T(0), // fitting of the atenuation, angular correction fSingleSectorFitter(0), // just for debugging // // Counters // fTotalTracks(0), // just for debugging - fAcceptedTracks(0), // just for debugging - fDebugCalPadRaw(0), // just for debugging - fDebugCalPadCorr(0), // just for debugging - fPrevIter(0) // the calibration object in its previous iteration (will not be owned by the new object, don't forget to delete it!) + fAcceptedTracks(0) // just for debugging { // // Constructor. // - G__SetCatchException(0); this->SetNameTitle(name, title); fCuts = cuts; - fPrevIter = prevIter; // // Fitter initialization // @@ -323,6 +380,13 @@ AliTPCcalibTracksGain::AliTPCcalibTracksGain(const char* name, const char* title fFitter1T = new TLinearFitter(45,"hyp44"); fFitter2T = new TLinearFitter(45,"hyp44"); // + fDFitter0M = new TLinearFitter(10,"hyp9"); + fDFitter1M = new TLinearFitter(10,"hyp9"); + fDFitter2M = new TLinearFitter(10,"hyp9"); + fDFitter0T = new TLinearFitter(10,"hyp9"); + fDFitter1T = new TLinearFitter(10,"hyp9"); + fDFitter2T = new TLinearFitter(10,"hyp9"); + // // fFitter0M->StoreData(kFALSE); fFitter1M->StoreData(kFALSE); @@ -331,64 +395,36 @@ AliTPCcalibTracksGain::AliTPCcalibTracksGain(const char* name, const char* title fFitter1T->StoreData(kFALSE); fFitter2T->StoreData(kFALSE); // + fDFitter0M->StoreData(kFALSE); + fDFitter1M->StoreData(kFALSE); + fDFitter2M->StoreData(kFALSE); + fDFitter0T->StoreData(kFALSE); + fDFitter1T->StoreData(kFALSE); + fDFitter2T->StoreData(kFALSE); // - // Add profile histograms -JUST for visualization - Not used for real calibration - // - // - fArrayQM=new TObjArray(73); // Qmax normalized - fArrayQT=new TObjArray(73); // Qtot normalized - fProfileArrayQM = new TObjArray(37); // Qmax normalized versus local X - fProfileArrayQT = new TObjArray(37); // Qtot normalized versus local X - fProfileArrayQM2D = new TObjArray(37); // Qmax normalized versus local X and phi - fProfileArrayQT2D = new TObjArray(37); // Qtot normalized versus local X and phi - char hname[1000]; - for (Int_t i=0; i<73; i++){ - sprintf(hname,"QM_%d",i); - fArrayQM->AddAt(new TH1F(hname,hname,200,0,1000),i); - sprintf(hname,"QT_%d",i); - fArrayQT->AddAt(new TH1F(hname,hname,200,0,1000),i); - } - - for (Int_t i=0; i<37;i++){ - sprintf(hname,"QMvsx_%d",i); - fProfileArrayQM->AddAt(new TProfile(hname,hname,50,89,250),i); - sprintf(hname,"QTvsx_%d",i); - fProfileArrayQT->AddAt(new TProfile(hname,hname,50,89,250),i); - sprintf(hname,"QM2D_%d",i); - fProfileArrayQM2D->AddAt(new TProfile2D(hname,hname,50,89,250,10,-0.15,0.15),i); - sprintf(hname,"QT2D_%d",i); - fProfileArrayQT2D->AddAt(new TProfile2D(hname,hname,50,89,250,10,-0.15,0.15),i); - } // // just for debugging -counters // fTotalTracks = 0; fAcceptedTracks = 0; - fDebugCalPadRaw = new AliTPCCalPad("DebugCalPadRaw", "All clusters simply added up before correction"); - fDebugCalPadCorr = new AliTPCCalPad("DebugCalPadCorr", "All clusters simply added up after correction"); - // this will be gone for the a new ROOT version > v5-17-05 - for (UInt_t i = 0; i < 36; i++) { - fNShortClusters[i] = 0; - fNMediumClusters[i] = 0; - fNLongClusters[i] = 0; - } - } +} AliTPCcalibTracksGain::~AliTPCcalibTracksGain() { // // Destructor. // - Info("Destructor",""); + Info("Destructor",":"); if (fSimpleFitter) delete fSimpleFitter; if (fSqrtFitter) delete fSqrtFitter; if (fLogFitter) delete fLogFitter; if (fSingleSectorFitter) delete fSingleSectorFitter; - if (fDebugCalPadRaw) delete fDebugCalPadRaw; - if (fDebugCalPadCorr) delete fDebugCalPadCorr; } + + + void AliTPCcalibTracksGain::Terminate(){ // // Evaluate fitters and close the debug stream. @@ -427,7 +463,14 @@ void AliTPCcalibTracksGain::Process(AliTPCseed* seed) { fFitter0T = new TLinearFitter(45,"hyp44"); fFitter1T = new TLinearFitter(45,"hyp44"); fFitter2T = new TLinearFitter(45,"hyp44"); - doinit=kTRUE; + // + fDFitter0M = new TLinearFitter(10,"hyp9"); + fDFitter1M = new TLinearFitter(10,"hyp9"); + fDFitter2M = new TLinearFitter(10,"hyp9"); + fDFitter0T = new TLinearFitter(10,"hyp9"); + fDFitter1T = new TLinearFitter(10,"hyp9"); + fDFitter2T = new TLinearFitter(10,"hyp9"); + doinit=kFALSE; } //} @@ -453,16 +496,13 @@ Long64_t AliTPCcalibTracksGain::Merge(TCollection *list) { if (!fSingleSectorFitter) fSingleSectorFitter = new AliTPCFitPad(8, "hyp7", ""); - // just for debugging - if (!fDebugCalPadRaw) fDebugCalPadRaw = new AliTPCCalPad("DebugCalPadRaw", "All clusters simply added up before correction"); - if (!fDebugCalPadCorr) fDebugCalPadCorr = new AliTPCCalPad("DebugCalPadCorr", "All clusters simply added up after correction"); TIterator* iter = list->MakeIterator(); AliTPCcalibTracksGain* cal = 0; while ((cal = (AliTPCcalibTracksGain*)iter->Next())) { if (!cal->InheritsFrom(AliTPCcalibTracksGain::Class())) { - Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName()); + //Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName()); return -1; } @@ -471,72 +511,73 @@ Long64_t AliTPCcalibTracksGain::Merge(TCollection *list) { return 0; } +Float_t AliTPCcalibTracksGain::GetGain(AliTPCclusterMI *cluster){ + // + // get gain + // + Float_t gainPad= 1; + AliTPCCalPad * gainMap = AliTPCcalibDB::Instance()->GetDedxGainFactor(); + if (gainMap) { + AliTPCCalROC * roc = gainMap->GetCalROC(cluster->GetDetector()); + gainPad = roc->GetValue(cluster->GetRow(), TMath::Nint(cluster->GetPad())); + } + return gainPad; +} + +Float_t AliTPCcalibTracksGain::GetMaxNorm(AliTPCclusterMI * cluster, Float_t ky, Float_t kz){ + // + // Get normalized amplituded if the gain map provided + // + AliTPCClusterParam * parcl = AliTPCcalibDB::Instance()->GetClusterParam(); + Float_t maxNorm = + parcl->QmaxCorrection(cluster->GetDetector(), cluster->GetRow(),cluster->GetPad(), + cluster->GetTimeBin(),ky,kz,0.5,0.5,1.6); + + return GetGain(cluster)*maxNorm; +} + + +Float_t AliTPCcalibTracksGain::GetQNorm(AliTPCclusterMI * cluster, Float_t ky, Float_t kz){ + // + // Get normalized amplituded if the gain map provided + // + AliTPCClusterParam * parcl = AliTPCcalibDB::Instance()->GetClusterParam(); + Float_t totNorm = parcl->QtotCorrection(cluster->GetDetector(), cluster->GetRow(),cluster->GetPad(), cluster->GetTimeBin(),ky,kz,0.5,0.5,cluster->GetQ(),2.5,1.6); + return GetGain(cluster)*totNorm; +} + + + void AliTPCcalibTracksGain::Add(AliTPCcalibTracksGain* cal) { // // Adds another AliTPCcalibTracksGain object to this object. // - fSimpleFitter->Add(cal->fSimpleFitter); - fSqrtFitter->Add(cal->fSqrtFitter); - fLogFitter->Add(cal->fLogFitter); - fSingleSectorFitter->Add(cal->fSingleSectorFitter); - // - // - // - fFitter0M->Add(cal->fFitter0M); - fFitter1M->Add(cal->fFitter1M); - fFitter2M->Add(cal->fFitter2M); - fFitter0T->Add(cal->fFitter0T); - fFitter1T->Add(cal->fFitter1T); - fFitter2T->Add(cal->fFitter2T); - // - // - // histograms - // - for (Int_t i=0; i<73; i++){ - TH1F *his,*hism; - his = (TH1F*)fArrayQM->At(i); - hism = (TH1F*)cal->fArrayQM->At(i); - if (his && hism) his->Add(hism); - his = (TH1F*)fArrayQT->At(i); - hism = (TH1F*)cal->fArrayQT->At(i); - if (his && hism) his->Add(hism); - } - // - // - for (Int_t i=0; i<37; i++){ - TProfile *his,*hism; - his = (TProfile*)fProfileArrayQM->At(i); - hism = (TProfile*)cal->fProfileArrayQM->At(i); - if (his && hism) his->Add(hism); - his = (TProfile*)fProfileArrayQT->At(i); - hism = (TProfile*)cal->fProfileArrayQT->At(i); - if (his && hism) his->Add(hism); - } - // - // - for (Int_t i=0; i<37; i++){ - TProfile2D *his,*hism; - his = (TProfile2D*)fProfileArrayQM2D->At(i); - hism = (TProfile2D*)cal->fProfileArrayQM2D->At(i); - if (his && hism) his->Add(hism); - his = (TProfile2D*)fProfileArrayQT2D->At(i); - hism = (TProfile2D*)cal->fProfileArrayQT2D->At(i); - if (his && hism) his->Add(hism); - } + fSimpleFitter->Add(cal->fSimpleFitter); + fSqrtFitter->Add(cal->fSqrtFitter); + fLogFitter->Add(cal->fLogFitter); + fSingleSectorFitter->Add(cal->fSingleSectorFitter); + // + // + // + if (cal->fFitter0M->GetNpoints()>0) fFitter0M->Add(cal->fFitter0M); + if (cal->fFitter1M->GetNpoints()>0) fFitter1M->Add(cal->fFitter1M); + if (cal->fFitter2M->GetNpoints()>0) fFitter2M->Add(cal->fFitter2M); + if (cal->fFitter0T->GetNpoints()>0) fFitter0T->Add(cal->fFitter0T); + if (cal->fFitter1T->GetNpoints()>0) fFitter1T->Add(cal->fFitter1T); + if (cal->fFitter2T->GetNpoints()>0) fFitter2T->Add(cal->fFitter2T); + // + if (cal->fDFitter0M->GetNpoints()>0) fDFitter0M->Add(cal->fDFitter0M); + if (cal->fDFitter1M->GetNpoints()>0) fDFitter1M->Add(cal->fDFitter1M); + if (cal->fDFitter2M->GetNpoints()>0) fDFitter2M->Add(cal->fDFitter2M); + if (cal->fDFitter0T->GetNpoints()>0) fDFitter0T->Add(cal->fDFitter0T); + if (cal->fDFitter1T->GetNpoints()>0) fDFitter1T->Add(cal->fDFitter1T); + if (cal->fDFitter2T->GetNpoints()>0) fDFitter2T->Add(cal->fDFitter2T); // - // this will be gone for the a new ROOT version > v5-17-05 - for (UInt_t iSegment = 0; iSegment < 36; iSegment++) { - fNShortClusters[iSegment] += cal->fNShortClusters[iSegment]; - fNMediumClusters[iSegment] += cal->fNMediumClusters[iSegment]; - fNLongClusters[iSegment] += cal->fNLongClusters[iSegment]; - } // just for debugging, remove me fTotalTracks += cal->fTotalTracks; fAcceptedTracks += cal->fAcceptedTracks; - fDebugCalPadRaw->Add(cal->fDebugCalPadRaw); - fDebugCalPadCorr->Add(cal->fDebugCalPadCorr); } @@ -547,66 +588,14 @@ void AliTPCcalibTracksGain::AddTrack(AliTPCseed* seed) { // DumpTrack(seed); - // - // simple histograming part - for (Int_t i=0; i<159; i++){ - AliTPCclusterMI* cluster = seed->GetClusterPointer(i); - if (cluster) AddCluster(cluster); - } } -void AliTPCcalibTracksGain::AddCluster(AliTPCclusterMI* cluster){ - // - // Adding cluster information to the simple histograms - // No correction, fittings are applied - // - Float_t kThreshold=5; - if (cluster->GetX()<=0) return; - if (cluster->GetQ()<=kThreshold) return; - // - - - Int_t sector = cluster->GetDetector(); - TH1F * his; - his = GetQT(sector); - if (his) his->Fill(cluster->GetQ()); - his = GetQT(-1); - if (his) his->Fill(cluster->GetQ()); - his = GetQM(sector); - if (his) his->Fill(cluster->GetMax()); - his = GetQM(-1); - if (his) his->Fill(cluster->GetMax()); - // - sector = sector%36; - TProfile *prof; - prof = GetProfileQT(sector); - if (prof) prof->Fill(cluster->GetX(),cluster->GetQ()); - prof = GetProfileQT(-1); - if (prof) prof->Fill(cluster->GetX(),cluster->GetQ()); - prof = GetProfileQM(sector); - if (prof) prof->Fill(cluster->GetX(),cluster->GetMax()); - prof = GetProfileQM(-1); - if (prof) prof->Fill(cluster->GetX(),cluster->GetMax()); - // - Float_t phi = cluster->GetY()/cluster->GetX(); - TProfile2D *prof2; - prof2 = GetProfileQT2D(sector); - if (prof2) prof2->Fill(cluster->GetX(),phi,cluster->GetQ()); - prof2 = GetProfileQT2D(-1); - if (prof2) prof2->Fill(cluster->GetX(),phi,cluster->GetQ()); - prof2 = GetProfileQM2D(sector); - if (prof2) prof2->Fill(cluster->GetX(),phi,cluster->GetMax()); - prof2 = GetProfileQM2D(-1); - if (prof2) prof2->Fill(cluster->GetX(),phi,cluster->GetMax()); - - // -} void AliTPCcalibTracksGain::AddCluster(AliTPCclusterMI* cluster, Float_t /*momenta*/, Float_t/* mdedx*/, Int_t padType, Float_t xcenter, TVectorD& dedxQ, TVectorD& /*dedxM*/, Float_t /*fraction*/, Float_t fraction2, Float_t dedge, - TVectorD& /*parY*/, TVectorD& /*parZ*/, TVectorD& meanPos) { + TVectorD& parY, TVectorD& parZ, TVectorD& meanPos) { // // Adds cluster to the appropriate fitter for later analysis. // The charge used for the fit is the maximum charge for this specific cluster or the @@ -618,14 +607,38 @@ void AliTPCcalibTracksGain::AddCluster(AliTPCclusterMI* cluster, Float_t /*momen // gets the square root of the charge, and the log fitter gets fgkM*(1+q/fgkM), where q is the original charge // and fgkM==25. // + Float_t kedge = 3; + Float_t kfraction = 0.7; + Int_t kinorm = 2; + + + // Where to put selection on threshold? + // Defined by the Q/dEdxT variable - see debug streamer: + // + // Debug stream variables: (Where tu cut ?) + // chain0->Draw("Cl.fQ/dedxQ.fElements[1]>>his(100,0,3)","fraction2<0.6&&dedge>3","",1000000); + // mean 1 sigma 0.25 + // chain0->Draw("Cl.fMax/dedxM.fElements[1]>>his(100,0,3)","fraction2<0.6&&dedge>3","",1000000) + // mean 1 sigma 0.25 + // chain0->Draw("Cl.fQ/dedxQ.fElements[2]>>his(100,0,3)","fraction2<0.7&&dedge>3","",1000000) + // mean 1 sigma 0.29 + // chain0->Draw("Cl.fMax/dedxM.fElements[2]>>his(100,0,3)","fraction2<0.7&&dedge>3","",1000000) + // mean 1 sigma 0.27 + // chain0->Draw("Cl.fQ/dedxQ.fElements[3]>>his(100,0,3)","fraction2<0.8&&dedge>3","",1000000) + // mean 1 sigma 0.32 + // + // chain0->Draw("Cl.fQ/dedxQ.fElements[4]>>his(100,0,3)","fraction2<0.9&&dedge>3","",1000000) + // mean 1 sigma 0.4 + + // Fraction choosen 0.7 if (!cluster) { Error("AddCluster", "Cluster not valid."); return; } - if (dedge < 3.) return; - if (fraction2 > 0.7) return; + if (dedge < kedge) return; + if (fraction2 > kfraction) return; //Int_t padType = GetPadType(cluster->GetX()); Double_t xx[7]; @@ -640,42 +653,35 @@ void AliTPCcalibTracksGain::AddCluster(AliTPCclusterMI* cluster, Float_t /*momen xx[4] = xx[0] * xx[1]; xx[5] = TMath::Abs(cluster->GetZ()) - TMath::Abs(meanPos[4]); xx[6] = xx[5] * xx[5]; - // - // Update profile histograms - // // // Update fitters // Int_t segment = cluster->GetDetector() % 36; - Double_t q = fgkUseTotalCharge ? ((Double_t)(cluster->GetQ())) : ((Double_t)(cluster->GetMax())); // note: no normalization to pad size! - - // just for debugging - Int_t row = 0; - Int_t pad = 0; - GetRowPad(cluster->GetX(), cluster->GetY(), row, pad); - fDebugCalPadRaw->GetCalROC(cluster->GetDetector())->SetValue(row, pad, q + fDebugCalPadRaw->GetCalROC(cluster->GetDetector())->GetValue(row, pad)); - + Double_t q = fgkUseTotalCharge ? + ((Double_t)(cluster->GetQ()/GetQNorm(cluster,parY[1], parZ[1]))) : ((Double_t)(cluster->GetMax()/GetMaxNorm(cluster,parY[1], parZ[1]))); + // correct charge by normalising to mean charge per track - q /= dedxQ[2]; + q /= dedxQ[kinorm]; // just for debugging - fDebugCalPadCorr->GetCalROC(cluster->GetDetector())->SetValue(row, pad, q + fDebugCalPadCorr->GetCalROC(cluster->GetDetector())->GetValue(row, pad)); Double_t sqrtQ = TMath::Sqrt(q); Double_t logQ = fgkM * TMath::Log(1 + q / fgkM); - fSimpleFitter->GetFitter(segment, padType)->AddPoint(xx, q); - fSqrtFitter->GetFitter(segment, padType)->AddPoint(xx, sqrtQ); - fLogFitter->GetFitter(segment, padType)->AddPoint(xx, logQ); - fSingleSectorFitter->GetFitter(0, padType)->AddPoint(xx, q); - - // this will be gone for the a new ROOT version > v5-17-05 - if (padType == kShortPads) - fNShortClusters[segment]++; - if (padType == kMediumPads) - fNMediumClusters[segment]++; - if (padType == kLongPads) - fNLongClusters[segment]++; + TLinearFitter * fitter =0; + // + fitter = fSimpleFitter->GetFitter(segment, padType); + fitter->AddPoint(xx, q); + // + fitter = fSqrtFitter->GetFitter(segment, padType); + fitter->AddPoint(xx, sqrtQ); + // + fitter = fLogFitter->GetFitter(segment, padType); + fitter->AddPoint(xx, logQ); + // + fitter=fSingleSectorFitter->GetFitter(0, padType); + fitter->AddPoint(xx, q); + } void AliTPCcalibTracksGain::Evaluate(Bool_t robust, Double_t frac) { @@ -696,6 +702,13 @@ void AliTPCcalibTracksGain::Evaluate(Bool_t robust, Double_t frac) { fFitter0T->Eval(); fFitter1T->Eval(); fFitter2T->Eval(); + // + fDFitter0M->Eval(); + fDFitter1M->Eval(); + fDFitter2M->Eval(); + fDFitter0T->Eval(); + fDFitter1T->Eval(); + fDFitter2T->Eval(); } AliTPCCalPad* AliTPCcalibTracksGain::CreateFitCalPad(UInt_t fitType, Bool_t undoTransformation, Bool_t normalizeToPadSize) { @@ -841,15 +854,24 @@ AliTPCCalROC* AliTPCcalibTracksGain::CreateCombinedCalROC(const AliTPCCalROC* ro return roc; } -void AliTPCcalibTracksGain::GetParameters(UInt_t segment, UInt_t padType, UInt_t fitType, TVectorD &fitParam) { +Bool_t AliTPCcalibTracksGain::GetParameters(UInt_t segment, UInt_t padType, UInt_t fitType, TVectorD &fitParam) { // // Puts the fit parameters for the specified segment (IROC & OROC), padType and fitType // into the fitParam TVectorD (which should contain 8 elements). // padType is one of kShortPads, kMediumPads, kLongPads. fitType is one of kSimpleFitter, kSqrtFitter, kLogFitter. // Note: The fitter has to be evaluated first! // - - GetFitter(segment, padType, fitType)->GetParameters(fitParam); + TLinearFitter * fitter = GetFitter(segment, padType, fitType); + if (fitter){ + fitter->Eval(); + fitter->GetParameters(fitParam); + return kTRUE; + }else{ + Error("AliTPCcalibTracksGain::GetParameters", + Form("Fitter%d_%d_%d not availble", segment, padType, fitType)); + return kFALSE; + } + return kFALSE; } void AliTPCcalibTracksGain::GetErrors(UInt_t segment, UInt_t padType, UInt_t fitType, TVectorD &fitError) { @@ -861,31 +883,9 @@ void AliTPCcalibTracksGain::GetErrors(UInt_t segment, UInt_t padType, UInt_t fit // GetFitter(segment, padType, fitType)->GetErrors(fitError); - fitError *= TMath::Sqrt(GetRedChi2(segment, padType, fitType)); + //fitError *= TMath::Sqrt(GetRedChi2(segment, padType, fitType)); } -Double_t AliTPCcalibTracksGain::GetRedChi2(UInt_t segment, UInt_t padType, UInt_t fitType) { - // - // Returns the reduced chi^2 value for the specified segment, padType and fitType. - // padType is one of kShortPads, kMediumPads, kLongPads. fitType is one of kSimpleFitter, kSqrtFitter, kLogFitter. - // Note: The fitter has to be evaluated first! - // - - // this will be gone for the a new ROOT version > v5-17-05 - Int_t lNClusters = 0; - switch (padType) { - case kShortPads: - lNClusters = fNShortClusters[segment]; - break; - case kMediumPads: - lNClusters = fNMediumClusters[segment]; - break; - case kLongPads: - lNClusters = fNLongClusters[segment]; - break; - } - return GetFitter(segment, padType, fitType)->GetChisquare()/(lNClusters - 8); -} void AliTPCcalibTracksGain::GetCovarianceMatrix(UInt_t segment, UInt_t padType, UInt_t fitType, TMatrixD& covMatrix) { // @@ -950,35 +950,6 @@ Int_t AliTPCcalibTracksGain::GetPadType(Double_t lx) { return -1; } -// ONLY FOR DEBUGGING PURPOSES - REMOVE ME WHEN NOT NEEDED ANYMORE -Bool_t AliTPCcalibTracksGain::GetRowPad(Double_t lx, Double_t ly, Int_t& row, Int_t& pad) { - // - // Calculate the row and pad number when the local coordinates are given. - // Returns kFALSE if the position is out of range, otherwise return kTRUE. - // WARNING: This function is preliminary and probably isn't very accurate!! - // - - Double_t irocLow = fgTPCparam->GetPadRowRadiiLow(0) - fgTPCparam->GetInnerPadPitchLength()/2; - //Double_t irocUp = fgTPCparam->GetPadRowRadiiLow(fgTPCparam->GetNRowLow()-1) + fgTPCparam->GetInnerPadPitchLength()/2; - Double_t orocLow1 = fgTPCparam->GetPadRowRadiiUp(0) - fgTPCparam->GetOuter1PadPitchLength()/2; - //Double_t orocUp1 = fgTPCparam->GetPadRowRadiiUp(fgTPCparam->GetNRowUp1()-1) + fgTPCparam->GetOuter1PadPitchLength()/2; - Double_t orocLow2 = fgTPCparam->GetPadRowRadiiUp(fgTPCparam->GetNRowUp1()) - fgTPCparam->GetOuter2PadPitchLength()/2; - //Double_t orocUp2 = fgTPCparam->GetPadRowRadiiUp(fgTPCparam->GetNRowUp()-1) + fgTPCparam->GetOuter2PadPitchLength()/2; - - if (GetPadType(lx) == 0) { - row = (Int_t)((lx - irocLow) / fgTPCparam->GetInnerPadPitchLength()); - pad = (Int_t)((ly + fgTPCparam->GetYInner(row)) / fgTPCparam->GetInnerPadPitchWidth()); - } else if (GetPadType(lx) == 1) { - row = (Int_t)((lx - orocLow1) / fgTPCparam->GetOuter1PadPitchLength()); - pad = (Int_t)((ly + fgTPCparam->GetYOuter(row)) / fgTPCparam->GetOuterPadPitchWidth()); - } else if (GetPadType(lx) == 2) { - row = fgTPCparam->GetNRowUp1() + (Int_t)((lx - orocLow2) / fgTPCparam->GetOuter2PadPitchLength()); - pad = (Int_t)((ly + fgTPCparam->GetYOuter(row)) / fgTPCparam->GetOuterPadPitchWidth()); - } - else return kFALSE; - return kTRUE; -} - void AliTPCcalibTracksGain::DumpTrack(AliTPCseed* track) { // // Dump track information to the debug stream @@ -1009,6 +980,11 @@ void AliTPCcalibTracksGain::DumpTrack(AliTPCseed* track) { TTreeSRedirector * cstream = GetDebugStreamer(); if (cstream){ (*cstream) << "Track" << + "run="<1 && count>1){ (*cstream) << "TrackG" << + "run="<GetNRowLow(); - xcenter = 108.47; - } - if (padType == 1) { - firstRow = fgTPCparam->GetNRowLow(); - lastRow = fgTPCparam->GetNRowLow() + fgTPCparam->GetNRowUp1(); - xcenter = 166.60; - } - if (padType == 2) { - firstRow = fgTPCparam->GetNRowLow() + fgTPCparam->GetNRowUp1(); - lastRow = fgTPCparam->GetNRowLow() + fgTPCparam->GetNRowUp(); - xcenter = 222.6; + Int_t firstRow = 0, lastRow = 0; + Int_t minRow = 100; + Float_t xcenter = 0; + const Float_t ktany = TMath::Tan(TMath::DegToRad() * 10); + const Float_t kedgey = 4.; + if (padType == 0) { + firstRow = 0; + lastRow = fgTPCparam->GetNRowLow(); + xcenter = 108.47; + } + if (padType == 1) { + firstRow = fgTPCparam->GetNRowLow(); + lastRow = fgTPCparam->GetNRowLow() + fgTPCparam->GetNRowUp1(); + xcenter = 166.60; } - minRow = (lastRow - firstRow) / 2; - // - // - Int_t nclusters = 0; - Int_t nclustersNE = 0; // number of not edge clusters - Int_t lastSector = -1; - Float_t amplitudeQ[100]; - Float_t amplitudeM[100]; + if (padType == 2) { + firstRow = fgTPCparam->GetNRowLow() + fgTPCparam->GetNRowUp1(); + lastRow = fgTPCparam->GetNRowLow() + fgTPCparam->GetNRowUp(); + xcenter = 222.6; + } + minRow = (lastRow - firstRow) / 2; + // + // + Int_t nclusters = 0; + Int_t nclustersNE = 0; // number of not edge clusters + Int_t lastSector = -1; + Float_t amplitudeQ[100]; + Float_t amplitudeM[100]; Int_t rowIn[100]; Int_t index[100]; // @@ -1105,8 +1088,8 @@ Bool_t AliTPCcalibTracksGain::GetDedx(AliTPCseed* track, Int_t padType, Int_t* / Int_t detector = cluster->GetDetector() ; if (lastSector == -1) lastSector = detector; if (lastSector != detector) continue; - amplitudeQ[nclusters] = cluster->GetQ(); - amplitudeM[nclusters] = cluster->GetMax(); + amplitudeQ[nclusters] = cluster->GetQ()/GetQNorm(cluster,parY[1], parZ[1]); + amplitudeM[nclusters] = cluster->GetMax()/GetMaxNorm(cluster,parY[1], parZ[1]); rowIn[nclusters] = iCluster; nclusters++; Double_t dx = cluster->GetX() - xcenter; @@ -1201,25 +1184,44 @@ Bool_t AliTPCcalibTracksGain::GetDedx(AliTPCseed* track, Int_t padType, Int_t* / Float_t fraction2 = Float_t(inonEdge) / Float_t(nclustersNE); AddCluster(cluster, momenta, mdedx, padType, xcenter, dedxQ, dedxM, fraction, fraction2, dedge, parY, parZ, meanPos); + + Double_t qNorm = GetQNorm(cluster,parY[1], parZ[1]); + Double_t mNorm = GetMaxNorm(cluster,parY[1], parZ[1]); + + Float_t gain = GetGain(cluster); if (cstream) (*cstream) << "dEdx" << - "Cl.=" << cluster << // cluster of interest - "P=" << momenta << // track momenta - "dedx=" << mdedx << // mean dedx - corrected for angle - "IPad=" << padType << // pad type 0..2 - "xc=" << xcenter << // x center of chamber - "dedxQ.=" << &dedxQ << // dedxQ - total charge - "dedxM.=" << &dedxM << // dedxM - maximal charge - "fraction=" << fraction << // fraction - order in statistic (0,1) - "fraction2=" << fraction2 << // fraction - order in statistic (0,1) - "dedge=" << dedge << // distance to the edge - "parY.=" << &parY << // line fit - "parZ.=" << &parZ << // line fit - "meanPos.=" << &meanPos << // mean position (dx, dx^2, y,y^2, z, z^2) - "\n"; + "run="<