]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TRD/AliTRDCalibTask.cxx
re-activate contrib code
[u/mrichter/AliRoot.git] / TRD / AliTRDCalibTask.cxx
index 7e4f6cc47a2c4552856f1051f8ec67b0e9c85bc8..15b67dd893f5578a408fe2e18544e30b231a41a4 100644 (file)
@@ -78,6 +78,9 @@ using namespace std;
 #include "AliLog.h"
 
 #include "AliTRDCalibChamberStatus.h"
+#include "AliTRDdEdxBaseUtils.h"
+#include "AliTRDdEdxCalibHistArray.h"
+#include "AliTRDdEdxCalibUtils.h"
 
 #include "AliTRDCalibTask.h"
 
@@ -85,82 +88,96 @@ using namespace std;
 ClassImp(AliTRDCalibTask)
 
 //________________________________________________________________________
-  AliTRDCalibTask::AliTRDCalibTask(const char *name) 
-    : AliAnalysisTaskSE(name), fESD(0),
-      fESDfriend(0),
-      fkEsdTrack(0),
-      fFriendTrack(0),
-      fCalibObject(0),
-      fTrdTrack(0),
-      fCl(0),
-      fListHist(0),
-      fTRDCalibraFillHisto(0),
-      fTRDChamberStatus(0),
-      fNEvents(0),
-      fNEventsInput(0),
-      fNbTRDTrack(0),
-      fNbTRDTrackOffline(0),
-      fNbTRDTrackStandalone(0),
-      fNbTPCTRDtrack(0),
-      fNbGoodTracks(0),
-      fNbTimeBin(0),
-      fNbTimeBinOffline(0),
-      fNbTimeBinStandalone(0),
-      fNbClusters(0),
-      fNbClustersOffline(0),
-      fNbClustersStandalone(0),
-      fNbTracklets(0),
-      fNbTrackletsOffline(0),
-      fNbTrackletsStandalone(0),
-      fAbsoluteGain(0),
-      fCH2dSum(0),
-      fPH2dSum(0),
-      fCH2dSM(0),
-      fPH2dSM(0),
-      fCH2dTest(0),
-      fPH2dTest(0),
-      fLinearVdriftTest(0),
-      fOnInstance(kTRUE),
-      fHisto2d(kTRUE),
-      fVector2d(kFALSE),
-      fVdriftLinear(kTRUE),
-      fExbAlt(kFALSE),
-      fNbTimeBins(0),
-      fSelectedTrigger(new TObjArray()),
-      fRejected(kTRUE),
-      fEsdTrackCuts(0),
-      fRequirePrimaryVertex(kFALSE),
-      fVtxTPC(kFALSE),
-      fVtxSPD(kFALSE),
-      fMinNbContributors(0),
-      fRangePrimaryVertexZ(9999999.0),
-      fMinNbTracks(9),
-      fMaxNbTracks(500),
-      fLow(0),
-      fHigh(30),
-      fFillZero(kFALSE),
-      fNormalizeNbOfCluster(kFALSE),
-      fRelativeScale(0.0),
-      fMaxCluster(100.0),
-      fNbMaxCluster(2),
-      fOfflineTracks(kFALSE),
-      fStandaloneTracks(kFALSE),
-      fFirstRunGain(-1),
-      fVersionGainUsed(-1),
-      fSubVersionGainUsed(-1),
-      fFirstRunGainLocal(-1),
-      fVersionGainLocalUsed(-1),
-      fSubVersionGainLocalUsed(-1),
-      fFirstRunVdrift(-1),
-      fVersionVdriftUsed(-1), 
-      fSubVersionVdriftUsed(-1),
-      fFirstRunExB(-1),
-      fVersionExBUsed(-1), 
-      fSubVersionExBUsed(-1),
-      fCalDetGain(0x0),
-      fMaxEvent(0),
-      fCounter(0),
-      fDebug(0)
+AliTRDCalibTask::AliTRDCalibTask(const char *name) 
+: AliAnalysisTaskSE(name), fESD(0),
+  fkEsdTrack(0),
+  fFriendTrack(0),
+  fCalibObject(0),
+  fTrdTrack(0),
+  fCl(0),
+  fListHist(0),
+  fTRDCalibraFillHisto(0),
+  fTRDChamberStatus(0),
+  fNEvents(0),
+  fNEventsInput(0),
+  fNbTRDTrack(0),
+  fNbTRDTrackOffline(0),
+  fNbTRDTrackStandalone(0),
+  fNbTPCTRDtrack(0),
+  fNbGoodTracks(0),
+  fNbTimeBin(0),
+  fNbTimeBinOffline(0),
+  fNbTimeBinStandalone(0),
+  fNbClusters(0),
+  fNbClustersOffline(0),
+  fNbClustersStandalone(0),
+  fNbTracklets(0),
+  fNbTrackletsOffline(0),
+  fNbTrackletsStandalone(0),
+  fAbsoluteGain(0),
+  fTOFbc(0),
+  fCH2dSum(0),
+  fPH2dSum(0),
+  fCH2dSM(0),
+  fPH2dSM(0),
+  fCH2dTest(0),
+  fPH2dTest(0),
+  fLinearVdriftTest(0),
+  fOnInstance(kTRUE),
+  fHisto2d(kTRUE),
+  fVector2d(kFALSE),
+  fVdriftLinear(kTRUE),
+  fExbAlt(kFALSE),
+  fDebugLevelTRDCalibraFillHisto(0),
+  fNbTimeBins(0),
+  fNumberBinCharge(100),
+  fRangeCharge(150),
+  fVdBindx(32), 
+  fVdBindy(70), 
+  fVdRangex(0.8),
+  fVdRangey(1.4),
+  fSelectTrigger(kTRUE),
+  fSelectedTrigger(new TObjArray()),
+  fRejected(kTRUE),
+  fEsdTrackCuts(0),
+  fRequirePrimaryVertex(kFALSE),
+  fVtxTPC(kFALSE),
+  fVtxSPD(kFALSE),
+  fMinNbContributors(0),
+  fRangePrimaryVertexZ(9999999.0),
+  fRejectPileUpWithSPD(kFALSE),
+  fRejectPileUpWithTOF(kFALSE),
+  fMinNbTracks(9),
+  fMaxNbTracks(999999999),
+  fCutWithVdriftCalib(kFALSE),
+  fMinNbTRDtracklets(0),
+  fMinTRDMomentum(0),
+  fScaleGainWithTPCSignal(kFALSE),
+  fLow(0),
+  fHigh(30),
+  fFillZero(kFALSE),
+  fNormalizeNbOfCluster(kFALSE),
+  fRelativeScale(0.0),
+  fMaxCluster(100.0),
+  fNbMaxCluster(2),
+  fOfflineTracks(kFALSE),
+  fStandaloneTracks(kFALSE),
+  fFirstRunGain(-1),
+  fVersionGainUsed(-1),
+  fSubVersionGainUsed(-1),
+  fFirstRunGainLocal(-1),
+  fVersionGainLocalUsed(-1),
+  fSubVersionGainLocalUsed(-1),
+  fFirstRunVdrift(-1),
+  fVersionVdriftUsed(-1), 
+  fSubVersionVdriftUsed(-1),
+  fFirstRunExB(-1),
+  fVersionExBUsed(-1), 
+  fSubVersionExBUsed(-1),
+  fCalDetGain(0x0),
+  fMaxEvent(0),
+  fCounter(0),
+  fPHQon(kTRUE)
 {
   //
   // Default constructor
@@ -174,6 +191,7 @@ ClassImp(AliTRDCalibTask)
   fNrphi[1] = 0;
   fNrphi[2] = 0;
 
+  
   // Define input and output slots here
   // Input slot #0 works with a TChain
   DefineInput(0, TChain::Class());
@@ -208,6 +226,7 @@ AliTRDCalibTask::~AliTRDCalibTask()
   if(fNbTrackletsOffline) delete fNbTrackletsOffline;
   if(fNbTrackletsStandalone) delete fNbTrackletsStandalone;
   if(fAbsoluteGain) delete fAbsoluteGain;
+  if(fTOFbc) delete fTOFbc;
   if(fCH2dSum) delete fCH2dSum;
   if(fPH2dSum) delete fPH2dSum;
   if(fCH2dSM) delete fCH2dSM;
@@ -215,6 +234,10 @@ AliTRDCalibTask::~AliTRDCalibTask()
   if(fCH2dTest) delete fCH2dTest;
   if(fPH2dTest) delete fPH2dTest;
   if(fLinearVdriftTest) delete fLinearVdriftTest;
+  if(IsPHQon()){
+    AliTRDdEdxCalibUtils::DeleteHistArray();
+  }
+
   if(fCalDetGain) delete fCalDetGain;
   
   if(fSelectedTrigger) {
@@ -259,6 +282,10 @@ void AliTRDCalibTask::UserCreateOutputObjects()
   // instance calibration 
   fTRDCalibraFillHisto = AliTRDCalibraFillHisto::Instance();
   if(fOnInstance) {
+    fTRDCalibraFillHisto->SetNumberBinCharge(fNumberBinCharge); // set number of bin of the charge distribution
+    fTRDCalibraFillHisto->SetCutWithVdriftCalib(fCutWithVdriftCalib); // cut vdrift
+    fTRDCalibraFillHisto->SetMinNbTRDtracklets(fMinNbTRDtracklets); // min number of TRD tracklets
+    fTRDCalibraFillHisto->SetMinTRDMomentum(fMinTRDMomentum); // min TRD momentum
     fTRDCalibraFillHisto->SetHisto2d(fHisto2d); // choose to use histograms
     fTRDCalibraFillHisto->SetVector2d(fVector2d); // choose to use vectors
     fTRDCalibraFillHisto->SetCH2dOn();  // choose to calibrate the gain
@@ -267,6 +294,7 @@ void AliTRDCalibTask::UserCreateOutputObjects()
     fTRDCalibraFillHisto->SetLinearFitterOn(fVdriftLinear); // Other possibility vdrift VDRIFT
     fTRDCalibraFillHisto->SetLinearFitterDebugOn(fVdriftLinear); // Other possibility vdrift
     fTRDCalibraFillHisto->SetExbAltFitOn(fExbAlt); // Alternative method for exb
+    fTRDCalibraFillHisto->SetScaleWithTPCSignal(fScaleGainWithTPCSignal); // Scale Gain with TPC signal
     for(Int_t k = 0; k < 3; k++){
       if(((fNz[k] != 10) && (fNrphi[k] != 10)) && ((fNz[k] != 100) && (fNrphi[k] != 100))) {
        fTRDCalibraFillHisto->SetNz(k,fNz[k]);                                    // Mode calibration
@@ -295,14 +323,23 @@ void AliTRDCalibTask::UserCreateOutputObjects()
     fTRDCalibraFillHisto->SetNumberClustersf(fHigh); // At least 11 clusters
     
     // For testing only
-    if(fDebug > 2) fTRDCalibraFillHisto->SetDebugLevel(1); //debug stuff
+    fTRDCalibraFillHisto->SetDebugLevel(fDebugLevelTRDCalibraFillHisto); //debug stuff
     
     if(fHisto2d) {  
       fListHist->Add(fTRDCalibraFillHisto->GetCH2d());
       fListHist->Add(fTRDCalibraFillHisto->GetPH2d()); 
       fListHist->Add(fTRDCalibraFillHisto->GetPRF2d());
     } 
-    if(fVdriftLinear) fListHist->Add((TObject *)fTRDCalibraFillHisto->GetVdriftLinearFit());
+    if(fVdriftLinear) {
+      AliTRDCalibraVdriftLinearFit *fvdl = fTRDCalibraFillHisto->GetVdriftLinearFit();
+      if(fvdl){
+       fvdl->SetNbBindx(fVdBindx);
+       fvdl->SetNbBindy(fVdBindy); 
+       fvdl->SetRangedx(fVdRangex); 
+       fvdl->SetRangedy(fVdRangey);  
+      }
+      fListHist->Add((TObject *)fTRDCalibraFillHisto->GetVdriftLinearFit());
+    }
     if(fVector2d) fListHist->Add((TObject *) fTRDCalibraFillHisto->GetCalibraVector()); //calibra vector
     if(fExbAlt) fListHist->Add((TObject *)fTRDCalibraFillHisto->GetExbAltFit());
   }
@@ -330,7 +367,24 @@ void AliTRDCalibTask::UserCreateOutputObjects()
   fAbsoluteGain->SetStats(0);
   fAbsoluteGain->Sumw2();
   fListHist->Add(fAbsoluteGain);
+
+  fTOFbc = new TH2F(Form("TOFbc_%s",(const char*)fName),"TOFbc", 200, -150., 50., nBinsPt, binLimPt);
+  fTOFbc->SetYTitle("Momentum");
+  fTOFbc->SetXTitle("TOF bc");
+  fTOFbc->SetZTitle("counts");
+  fTOFbc->SetStats(0);
+  fTOFbc->Sumw2();
+  fListHist->Add(fTOFbc);
   
+  if(IsPHQon()){
+    printf("\n        AliTRDCalibTask PHQ is on!!     \n\n");
+    AliTRDdEdxBaseUtils::PrintControl();
+    AliTRDdEdxCalibUtils::IniHistArray(fListHist, kTRUE);
+  }
+  else{
+    printf("\n        AliTRDCalibTask PHQ is off!!     \n\n");
+  }
+
   /////////////////////////////////////////
   // First debug level
   ///////////////////////////////////////
@@ -552,13 +606,13 @@ void AliTRDCalibTask::UserExec(Option_t *)
   //cout << "event = " << fCounter << endl;
   
   //printf("Counter %d\n",fCounter);
-  
+
   ///////////////////
   // Check trigger
   ///////////////////
   Bool_t pass = kTRUE;
 
-  if (strstr(type,"p-p")) {
+  if (fSelectTrigger) {
    
     //printf("Will check the triggers\n");
 
@@ -620,7 +674,19 @@ void AliTRDCalibTask::UserExec(Option_t *)
   }
   
   //printf("Primary vertex passed\n");
+
+  //////////////////////////////////
+  // Reject pile-up with SPD
+  //////////////////////////////////
   
+  if(fRejectPileUpWithSPD) {
+    if(fInputEvent->IsPileupFromSPD(3, 0.8, 3., 2., 5)){
+      //printf("test\n");
+      PostData(1, fListHist);
+      return;
+    }
+  }
+
   //////////////////////////////////////
   // Requirement on number of good tracks
   //////////////////////////////////////
@@ -638,8 +704,10 @@ void AliTRDCalibTask::UserExec(Option_t *)
     //printf("centrality %f, centralityb %f\n",centrality,centralityb);
   }
   
-  if (strstr(type,"Pb-Pb")) {
-    //printf("Will check the number of good tracks\n");
+  //printf("Beam type %s\n",(const char*)type);
+  if (!strstr(type,"p-p")) {
+    //if (strstr(type,"A-A")) {
+    //printf("Will check the number of good tracks %f %f %f\n",nGoodParticles,fMinNbTracks,fMaxNbTracks);
     if((nGoodParticles < fMinNbTracks) || (nGoodParticles > fMaxNbTracks)) {
       PostData(1, fListHist);
       return;
@@ -658,7 +726,7 @@ void AliTRDCalibTask::UserExec(Option_t *)
   Int_t nbtrackTPC = 0;
   
 
-  
+
   if (nbTracks <= 0.0) {
     
     if(fDebug > 1) {
@@ -669,8 +737,8 @@ void AliTRDCalibTask::UserExec(Option_t *)
     PostData(1, fListHist);
     return;
   }
-  
-  
+
+
   fESDfriend = dynamic_cast<AliESDfriend*> (fESD->FindListObject("AliESDfriend"));
   if(!fESDfriend){
     AliError("fESDfriend not available");
@@ -682,14 +750,15 @@ void AliTRDCalibTask::UserExec(Option_t *)
     PostData(1, fListHist);
     return;
   }
-  
+
   //printf("has friends\n");
 
   /////////////////////////////////////
   // Loop on AliESDtrack
   ////////////////////////////////////
-  //printf("Nb of tracks %f\n",nbTracks);      
-  for(int itrk=0; itrk < nbTracks; ++itrk){
+  //printf("Nb of tracks %f\n",nbTracks); 
+  Int_t nbTracksfriends = fESDfriend->GetNumberOfTracks();     
+  for(int itrk=0; itrk < nbTracksfriends; ++itrk){
     
     // Get ESD track
     fkEsdTrack = fESD->GetTrack(itrk);
@@ -697,7 +766,10 @@ void AliTRDCalibTask::UserExec(Option_t *)
     ULong_t status = fkEsdTrack->GetStatus(); 
     if(status&(AliESDtrack::kTPCout)) ++nbtrackTPC;
     
+    // Fix suggested by Alex
     fFriendTrack = fESDfriend->GetTrack(itrk);
+    //printf("itrk %d\n",itrk);
+    //fFriendTrack = (fESDfriend->GetNumberOfTracks()>itrk)?fESDfriend->GetTrack(itrk):NULL;
     if(!fFriendTrack)  {
       //printf("No friend track %d\n",itrk);
       continue;
@@ -749,6 +821,13 @@ void AliTRDCalibTask::UserExec(Option_t *)
       continue;
     }
     
+    Int_t nbcrossing = fkEsdTrack->GetTOFBunchCrossing();
+    fTOFbc->Fill(nbcrossing,fkEsdTrack->Pt());
+    // TOF pile-up rejection is asked
+    if(fRejectPileUpWithTOF) {
+      if(TMath::Abs(nbcrossing)>0.5) continue;
+    }
+
     // First Absolute gain calibration
     Int_t trdNTracklets = (Int_t) fkEsdTrack->GetTRDntracklets();
     Int_t trdNTrackletsPID = (Int_t) fkEsdTrack->GetTRDntrackletsPID(); 
@@ -770,12 +849,19 @@ void AliTRDCalibTask::UserExec(Option_t *)
 
     if(good && fOnInstance) {
       //cout << "good" << endl;
-      fTRDCalibraFillHisto->UpdateHistogramsV1(fTrdTrack);
+      fTRDCalibraFillHisto->UpdateHistogramsV1(fTrdTrack,fkEsdTrack);
       //printf("Fill fTRDCalibraFillHisto\n");
     }
-      
-      
-         
+
+    if(IsPHQon()){
+      const Double_t mag = AliTRDdEdxBaseUtils::IsExBOn() ? fESD->GetMagneticField() : -1;
+      const Int_t charge = AliTRDdEdxBaseUtils::IsExBOn() ? fkEsdTrack->Charge() : -1;
+      const Double_t toTPCscale = AliTRDdEdxCalibUtils::GetCalibTPCscale(fkEsdTrack->GetTPCncls(), fkEsdTrack->GetTPCsignal());
+      if(toTPCscale>0){
+        AliTRDdEdxCalibUtils::FillHist(fTrdTrack, 0, mag, charge, toTPCscale);
+      }
+    }
+
     //////////////////////////////////
     // Debug 
     ////////////////////////////////
@@ -787,7 +873,8 @@ void AliTRDCalibTask::UserExec(Option_t *)
       Int_t nbtracklets = 0;
       
       // Check some stuff
-      Bool_t standalonetracklet = kFALSE;  
+      // Coverity
+      //Bool_t standalonetracklet = kFALSE;  
       const AliTRDseedV1 *tracklet = 0x0;
       //////////////////////////////////////
       // Loop tracklets
@@ -796,7 +883,7 @@ void AliTRDCalibTask::UserExec(Option_t *)
       Double_t phtb[AliTRDseedV1::kNtb];
       memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
       Double_t sum = 0.0;
-      Float_t normalisation = 6.67;
+      Float_t normalisation = 1.13;
       Int_t detector = 0;
       Int_t sector = 0;
       for(Int_t itr = 0; itr < 6; ++itr){
@@ -804,8 +891,9 @@ void AliTRDCalibTask::UserExec(Option_t *)
        if(!(tracklet = fTrdTrack->GetTracklet(itr))) continue;
        if(!tracklet->IsOK()) continue;
        ++nbtracklets;
-       standalonetracklet = kFALSE; 
-       if(tracklet->IsStandAlone()) standalonetracklet = kTRUE;
+        // Coverity
+       //standalonetracklet = kFALSE; 
+       //if(tracklet->IsStandAlone()) standalonetracklet = kTRUE;
        
        nbclusters = 0;
        memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
@@ -828,7 +916,8 @@ void AliTRDCalibTask::UserExec(Option_t *)
          if(!(fCl = tracklet->GetClusters(ic))) continue;
          ++nbclusters;
          time = fCl->GetPadTime();
-         ch =  tracklet->GetdQdl(ic);
+         //ch =  tracklet->GetdQdl(ic);
+         ch =  tracklet->GetQperTB(ic);
          qcl = TMath::Abs(fCl->GetQ());
          detector = fCl->GetDetector();          
          // Add the charge if shared cluster
@@ -1313,7 +1402,7 @@ Bool_t AliTRDCalibTask::SetVersionSubversion(){
   // Load Chamber Gain factors into the Tender supply
   //
   
-  printf("SetVersionSubversion\n");
+  //printf("SetVersionSubversion\n");
 
   //find previous entry from the UserInfo
   TTree *tree=((TChain*)GetInputData(0))->GetTree();
@@ -1356,12 +1445,13 @@ Bool_t AliTRDCalibTask::SetVersionSubversion(){
       fFirstRunGainLocal = id->GetFirstRun();
       fVersionGainLocalUsed = id->GetVersion();
       fSubVersionGainLocalUsed = id->GetSubVersion();
-    } else if(os->GetString().Contains("TRD/Calib/ChamberExB")){
+    } else if((os->GetString().Contains("TRD/Calib/ChamberExB")) && (!os->GetString().Contains("TRD/Calib/ChamberExBAlt"))){
       // Get Old drift velocity calibration
       AliCDBId *id=AliCDBId::MakeFromString(os->GetString());
       fFirstRunExB = id->GetFirstRun();
       fVersionExBUsed = id->GetVersion();
       fSubVersionExBUsed = id->GetSubVersion();
+      //printf("Version %d and subversion %d\n",fVersionExBUsed,fSubVersionExBUsed);
     }
   }
 
@@ -1423,7 +1513,8 @@ Bool_t AliTRDCalibTask::FindP1TrackPHtrackletV1Test(const AliTRDseedV1 *tracklet
   ////////////
   // results of the linear fit
   Double_t dydt                       = 0.0;                                // dydt tracklet after straight line fit
-  Double_t errorpar                   = 0.0;                                // error after straight line fit on dy/dt
+  // Coverity
+  //Double_t errorpar                   = 0.0;                                // error after straight line fit on dy/dt
   Double_t pointError                 = 0.0;                                // error after straight line fit 
   // pad row problemes: avoid tracklet that cross pad rows, tilting angle in the constant
   Int_t    crossrow                   = 0;                                  // if it crosses a pad row
@@ -1497,7 +1588,8 @@ Bool_t AliTRDCalibTask::FindP1TrackPHtrackletV1Test(const AliTRDseedV1 *tracklet
   linearFitterTracklet.Eval();
   linearFitterTracklet.GetParameters(pars);
   pointError  =  TMath::Sqrt(linearFitterTracklet.GetChisquare()/(nbli-2));
-  errorpar    =  linearFitterTracklet.GetParError(1)*pointError;
+  // Coverity
+  //errorpar    =  linearFitterTracklet.GetParError(1)*pointError;
   dydt        = pars[1]; 
   //printf("chis %f, nbli %d, pointError %f, parError %f, errorpar %f\n",linearFitterTracklet->GetChisquare(),nbli,pointError,linearFitterTracklet->GetParError(1),errorpar);
   linearFitterTracklet.ClearPoints();