// 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 <TPDGCode.h>
#include <TStyle.h>
#include "TSystem.h"
#include "AliTPCClusterParam.h"
#include "AliTrackPointArray.h"
#include "TCint.h"
-#include "AliTPCcalibTracksGain.h"
#include <TH1.h>
#include <TH3F.h>
#include <TLinearFitter.h>
#include <TFile.h>
#include <TCollection.h>
#include <TIterator.h>
+#include <TProfile.h>
+#include <TProfile2D.h>
+#include <TProof.h>
+#include <TStatToolkit.h>
//
// AliRoot includes
#include "AliTPCCalROC.h"
#include "AliTPCCalPad.h"
#include "AliTPCClusterParam.h"
+#include "AliTPCcalibDB.h"
//
#include "AliTracker.h"
#include "AliESD.h"
#include "AliTPCclusterMI.h"
#include "AliTPCcalibTracksCuts.h"
#include "AliTPCFitPad.h"
+#include "TStatToolkit.h"
+#include "TString.h"
+#include "TCut.h"
-// REMOVE ALL OF THIS
+//
#include <TTree.h>
#include "AliESDEvent.h"
AliTPCcalibTracksGain::AliTPCcalibTracksGain() :
AliTPCcalibBase(),
- fDebugStream(0), //! debug stream for debugging
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
{
//
AliTPCcalibTracksGain::AliTPCcalibTracksGain(const AliTPCcalibTracksGain& obj) :
AliTPCcalibBase(obj),
- fDebugStream(0), //! debug stream for debugging
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
//
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
{
//
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(),
- fDebugStream(0), //! debug stream for debugging
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
//
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
//
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");
//
- // Add profile histograms -JUST for visualization - Not used for real calibration
//
+ fFitter0M->StoreData(kFALSE);
+ fFitter1M->StoreData(kFALSE);
+ fFitter2M->StoreData(kFALSE);
+ fFitter0T->StoreData(kFALSE);
+ fFitter1T->StoreData(kFALSE);
+ fFitter2T->StoreData(kFALSE);
+ //
+ fDFitter0M->StoreData(kFALSE);
+ fDFitter1M->StoreData(kFALSE);
+ fDFitter2M->StoreData(kFALSE);
+ fDFitter0T->StoreData(kFALSE);
+ fDFitter1T->StoreData(kFALSE);
+ fDFitter2T->StoreData(kFALSE);
//
- 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 (fDebugStream) {
- //fDebugStream->GetFile()->Close();
- printf("Deleting debug stream object\n");
- delete fDebugStream;
- }
- if (fDebugCalPadRaw) delete fDebugCalPadRaw;
- if (fDebugCalPadCorr) delete fDebugCalPadCorr;
-}
void AliTPCcalibTracksGain::Terminate(){
//
//
Evaluate();
- if (fDebugStream) {
- delete fDebugStream;
- fDebugStream = 0;
- }
+ AliTPCcalibBase::Terminate();
}
// is added.
//
+
fTotalTracks++;
if (!fCuts->AcceptTrack(seed)) return;
+ //
+ // reinint on proof
+ // if (gProof){
+ static Bool_t doinit= kTRUE;
+ if (doinit){
+ fSimpleFitter = new AliTPCFitPad(8, "hyp7", "");
+ fSqrtFitter = new AliTPCFitPad(8, "hyp7", "");
+ fLogFitter = new AliTPCFitPad(8, "hyp7", "");
+ fSingleSectorFitter = new AliTPCFitPad(8, "hyp7", "");
+ //
+ fFitter0M = new TLinearFitter(45,"hyp44");
+ fFitter1M = new TLinearFitter(45,"hyp44");
+ fFitter2M = new TLinearFitter(45,"hyp44");
+ fFitter0T = new TLinearFitter(45,"hyp44");
+ 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");
+ doinit=kFALSE;
+ }
+ //}
+
fAcceptedTracks++;
AddTrack(seed);
}
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;
}
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++) {
// just for debugging, remove me
fTotalTracks += cal->fTotalTracks;
fAcceptedTracks += cal->fAcceptedTracks;
- fDebugCalPadRaw->Add(cal->fDebugCalPadRaw);
- fDebugCalPadCorr->Add(cal->fDebugCalPadCorr);
}
// See AddCluster(...) for more detail.
//
- if (!fDebugStream) fDebugStream = new TTreeSRedirector(fgkDebugStreamFileName);
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
// 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];
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);
+ 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);
// this will be gone for the a new ROOT version > v5-17-05
if (padType == kShortPads)
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) {
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) {
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
AddTracklet(sector[ipad],ipad, dedxQ[ipad], dedxM[ipad], parY[ipad], parZ[ipad], meanPos[ipad] );
if (isOK) count++;
}
-
- (*fDebugStream) << "Track" <<
- "Track.=" << track << // track information
- "\n";
- //
- //
- //
- if (count>1){
- (*fDebugStream) << "TrackG" <<
+
+ TTreeSRedirector * cstream = GetDebugStreamer();
+ if (cstream){
+ (*cstream) << "Track" <<
+ "run="<<fRun<< // run number
+ "event="<<fEvent<< // event number
+ "time="<<fTime<< // time stamp of event
+ "trigger="<<fTrigger<< // trigger
+ "mag="<<fMagF<< // magnetic field
"Track.=" << track << // track information
- "ncount="<<count<<
- // info for pad type 0
- "sector0="<<sector[0]<<
- "npoints0="<<npoints[0]<<
- "dedxM0.="<<&dedxM[0]<<
- "dedxQ0.="<<&dedxQ[0]<<
- "parY0.="<<&parY[0]<<
- "parZ0.="<<&parZ[0]<<
- "meanPos0.="<<&meanPos[0]<<
- //
- // info for pad type 1
- "sector1="<<sector[1]<<
- "npoints1="<<npoints[1]<<
- "dedxM1.="<<&dedxM[1]<<
- "dedxQ1.="<<&dedxQ[1]<<
- "parY1.="<<&parY[1]<<
- "parZ1.="<<&parZ[1]<<
- "meanPos1.="<<&meanPos[1]<<
- //
- // info for pad type 2
- "sector2="<<sector[2]<<
- "npoints2="<<npoints[2]<<
- "dedxM2.="<<&dedxM[2]<<
- "dedxQ2.="<<&dedxQ[2]<<
- "parY2.="<<&parY[2]<<
- "parZ2.="<<&parZ[2]<<
- "meanPos2.="<<&meanPos[2]<<
- //
"\n";
+ //
+ //
+ //
+ if ( GetStreamLevel()>1 && count>1){
+ (*cstream) << "TrackG" <<
+ "run="<<fRun<< // run number
+ "event="<<fEvent<< // event number
+ "time="<<fTime<< // time stamp of event
+ "trigger="<<fTrigger<< // trigger
+ "mag="<<fMagF<< // magnetic field
+ "Track.=" << track << // track information
+ "ncount="<<count<<
+ // info for pad type 0
+ "sector0="<<sector[0]<<
+ "npoints0="<<npoints[0]<<
+ "dedxM0.="<<&dedxM[0]<<
+ "dedxQ0.="<<&dedxQ[0]<<
+ "parY0.="<<&parY[0]<<
+ "parZ0.="<<&parZ[0]<<
+ "meanPos0.="<<&meanPos[0]<<
+ //
+ // info for pad type 1
+ "sector1="<<sector[1]<<
+ "npoints1="<<npoints[1]<<
+ "dedxM1.="<<&dedxM[1]<<
+ "dedxQ1.="<<&dedxQ[1]<<
+ "parY1.="<<&parY[1]<<
+ "parZ1.="<<&parZ[1]<<
+ "meanPos1.="<<&meanPos[1]<<
+ //
+ // info for pad type 2
+ "sector2="<<sector[2]<<
+ "npoints2="<<npoints[2]<<
+ "dedxM2.="<<&dedxM[2]<<
+ "dedxQ2.="<<&dedxQ[2]<<
+ "parY2.="<<&parY[2]<<
+ "parZ2.="<<&parZ[2]<<
+ "meanPos2.="<<&meanPos[2]<<
+ //
+ "\n";
+ }
}
-
}
Bool_t AliTPCcalibTracksGain::GetDedx(AliTPCseed* track, Int_t padType, Int_t* /*rows*/,
static TLinearFitter fitY(2, "pol1");
static TLinearFitter fitZ(2, "pol1");
+ fitY.StoreData(kFALSE);
+ fitZ.StoreData(kFALSE);
//
//
- 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;
- }
- 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];
//
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;
dedxQ[i] /= ndedx[i];
dedxM[i] /= ndedx[i];
}
-
+ TTreeSRedirector * cstream = GetDebugStreamer();
inonEdge = 0;
Float_t momenta = track->GetP();
Float_t mdedx = track->GetdEdx();
Float_t fraction2 = Float_t(inonEdge) / Float_t(nclustersNE);
AddCluster(cluster, momenta, mdedx, padType, xcenter, dedxQ, dedxM, fraction, fraction2, dedge, parY, parZ, meanPos);
-
- (*fDebugStream) << "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";
+
+ 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" <<
+ "run="<<fRun<< // run number
+ "event="<<fEvent<< // event number
+ "time="<<fTime<< // time stamp of event
+ "trigger="<<fTrigger<< // trigger
+ "mag="<<fMagF<< // magnetic field
+
+ "Cl.=" << cluster << // cluster of interest
+ "gain="<<gain<< // gain at cluster position
+ "mNorm="<<mNorm<< // Q max normalization
+ "qNorm="<<qNorm<< // Q tot normalization
+ "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";
}
- (*fDebugStream) << "dEdxT" <<
+ if (cstream) (*cstream) << "dEdxT" <<
+ "run="<<fRun<< // run number
+ "event="<<fEvent<< // event number
+ "time="<<fTime<< // time stamp of event
+ "trigger="<<fTrigger<< // trigger
+ "mag="<<fMagF<< // magnetic field
"P=" << momenta << // track momenta
"npoints="<<inonEdge<< // number of points
"sector="<<lastSector<< // sector number
if (padType==1) fitterT = fFitter1T;
if (padType==2) fitterT = fFitter2T;
fitterT->AddPoint(xxx,dedxQ[1]);
+ //
+ TLinearFitter *dfitterM = fDFitter0M;
+ if (padType==1) dfitterM=fDFitter1M;
+ if (padType==2) dfitterM=fDFitter2M;
+ dfitterM->AddPoint(xxx,dedxM[1]);
+ //
+ TLinearFitter *dfitterT = fDFitter0T;
+ if (padType==1) dfitterT = fDFitter1T;
+ if (padType==2) dfitterT = fDFitter2T;
+ dfitterT->AddPoint(xxx,dedxQ[1]);
}
// void SetQnorm(Int_t ipad, Int_t itype, TVectorD * norm);
TVectorD vec;
+
//
- fFitter0T->GetParameters(vec);
+ fDFitter0T->Eval();
+ fDFitter1T->Eval();
+ fDFitter2T->Eval();
+ fDFitter0M->Eval();
+ fDFitter1M->Eval();
+ fDFitter2M->Eval();
+ fDFitter0T->GetParameters(vec);
clparam->SetQnorm(0,0,&vec);
- fFitter1T->GetParameters(vec);
+ fDFitter1T->GetParameters(vec);
clparam->SetQnorm(1,0,&vec);
- fFitter2T->GetParameters(vec);
+ fDFitter2T->GetParameters(vec);
clparam->SetQnorm(2,0,&vec);
//
- fFitter0M->GetParameters(vec);
+ fDFitter0M->GetParameters(vec);
clparam->SetQnorm(0,1,&vec);
- fFitter1M->GetParameters(vec);
+ fDFitter1M->GetParameters(vec);
clparam->SetQnorm(1,1,&vec);
- fFitter2M->GetParameters(vec);
+ fDFitter2M->GetParameters(vec);
clparam->SetQnorm(2,1,&vec);
//
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff