Update of calibration code (Raphaelle)
authorcblume <cblume@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 2 Jun 2010 13:57:42 +0000 (13:57 +0000)
committercblume <cblume@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 2 Jun 2010 13:57:42 +0000 (13:57 +0000)
TRD/AliTRDCalibTask.cxx
TRD/AliTRDCalibTask.h
TRD/AliTRDCalibraFillHisto.cxx
TRD/AliTRDCalibraFillHisto.h
TRD/AliTRDPreprocessorOffline.cxx
TRD/AliTRDPreprocessorOffline.h
TRD/Macros/AddTaskTRDCalib.C

index a10cbc3..dafebb0 100644 (file)
@@ -25,8 +25,8 @@
 //                            
 //////////////////////////////////////////////////////////////////////////////////////
 
-
-
+#include <iostream>
+using namespace std;
 #include "Riostream.h"
 #include "TChain.h"
 #include "TTree.h"
@@ -77,7 +77,7 @@ ClassImp(AliTRDCalibTask)
 
 //________________________________________________________________________
   AliTRDCalibTask::AliTRDCalibTask(const char *name) 
-    : AliAnalysisTask(name,""), fESD(0),
+    : AliAnalysisTaskSE(name), fESD(0),
       fESDfriend(0),
       fkEsdTrack(0),
       fFriendTrack(0),
@@ -148,7 +148,7 @@ ClassImp(AliTRDCalibTask)
   DefineInput(0, TChain::Class());
         
   // Output slot #0 writes into a TList container
-  DefineOutput(0, TList::Class());
+  DefineOutput(1, TList::Class());
   
    
 }
@@ -190,37 +190,48 @@ AliTRDCalibTask::~AliTRDCalibTask()
   }
   
 }
+
+/*
 //________________________________________________________________________
 void AliTRDCalibTask::ConnectInputData(Option_t *) 
 {
   // Connect ESD or AOD here
-  // Called once
+  // Called once per event
   
-  TTree* tree = dynamic_cast<TTree*> (GetInputData(0)); //pointer wird "umgecastet" auf anderen Variablentyp
-  if (!tree) {
+  cout << "AliTRDCalibTask::ConnectInputData() IN" << endl;
+
+
+  //  TTree* tree = dynamic_cast<TTree*> (GetInputData(0)); //pointer wird "umgecastet" auf anderen Variablentyp
+  //  if (!tree) {
     //Printf("ERROR: Could not read chain from input slot 0");
-  } else {
-    
-    AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
-    
-    if (!esdH) {
-      //Printf("ERROR: Could not get ESDInputHandler");
-    } else {
-      fESD = esdH->GetEvent();
-      esdH->SetActiveBranches("ESDfriend*");
-      if ((esdH->GetTree())->GetBranch("ESDfriend.")) fESDfriend = esdH->GetESDfriend();
-      //else printf("No friend ESD\n");
-      //Printf("*** CONNECTED NEW EVENT ****");
-    }
+  //  } else {
     
+  AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
+  
+  if (!esdH) {
+    //Printf("ERROR: Could not get ESDInputHandler");
+  } else {
+    fESD = esdH->GetEvent();
+    //    esdH->SetActiveBranches("ESDfriend*");
+    if ((esdH->GetTree())->GetBranch("ESDfriend.")) fESDfriend = esdH->GetESDfriend();
+    //else printf("No friend ESD\n");
+    //Printf("*** CONNECTED NEW EVENT ****");
   }
+    
+
+    //  }
+  //cout << "AliTRDCalibTask::ConnectInputData() OUT" << endl;
+
 }
+*/
+
 //________________________________________________________________________
-void AliTRDCalibTask::CreateOutputObjects() 
+void AliTRDCalibTask::UserCreateOutputObjects() 
 {
   //
   // Create the histos
   //
+  //cout << "AliTRDCalibTask::CreateOutputObjects() IN" << endl;
 
   // Number of time bins
   if(fNbTimeBins==0) {
@@ -398,119 +409,134 @@ void AliTRDCalibTask::CreateOutputObjects()
     fListHist->Add(fNbTrackletsStandalone);
     
   }
-  
- }
- //________________________________________________________________________
- void AliTRDCalibTask::Exec(Option_t *) 
- {
-   //
-   // Filling of the histos
-   //
-
-   AliLog::SetGlobalLogLevel(AliLog::kError);
-   
-   if (!fESD) {
-     //Printf("ERROR: fESD not available");
-     return;
-   }
+  //cout << "AliTRDCalibTask::UserCreateOutputObjects() OUT" << endl;
+}
 
-   if (!fESDfriend) {
-     fESDfriend = (AliESDfriend*)(fESD->FindListObject("AliESDfriend"));
-     if(!fESDfriend) Printf("ERROR: fESDfriend not available");
-     //Printf("ERROR: fESDfriend not available");
-     return;
-   }
 
-   //printf("Counter %d\n",fCounter);
-
-   fCounter++;
-   if((fMaxEvent != 0) && (fMaxEvent < fCounter)) return;
-
-   //printf("Counter %d\n",fCounter);
-
-   ///////////////////
-   // Check trigger
-   ///////////////////
-   Bool_t pass = kTRUE;
-   Int_t numberOfTriggerSelected = fSelectedTrigger->GetEntriesFast();
-   //printf("numberofTriggerSelected %d\n",numberOfTriggerSelected);
-   if(fRejected) {
-     pass = kTRUE;
-     for(Int_t k = 0; k < numberOfTriggerSelected; k++){
-       const TObjString *const obString=(TObjString*)fSelectedTrigger->At(k);
-       const TString tString=obString->GetString();
-       if(fESD->IsTriggerClassFired((const char*)tString)) {
-        pass = kFALSE;
-       }
-     }
-   }
-   else {
-     pass = kFALSE;
-     for(Int_t k = 0; k < numberOfTriggerSelected; k++){
-       const TObjString *const obString=(TObjString*)fSelectedTrigger->At(k);
-       const TString tString=obString->GetString();
-       if(fESD->IsTriggerClassFired((const char*)tString)) {
-        pass = kTRUE;
-       }
-     }
-   }
-   if(!pass) return;   
-   //printf("Class Fired %s\n",(const char*)fESD->GetFiredTriggerClasses());
-   //printf("Trigger passed\n");
+//________________________________________________________________________
+void AliTRDCalibTask::UserExec(Option_t *) 
+{
+  //
+  // Filling of the histos
+  //
+  //cout << "AliTRDCalibTask::Exec() IN" << endl;
   
-   ///////////////////////////////
-   // Require a primary vertex
-   //////////////////////////////
-   if(fRequirePrimaryVertex) {
-     const AliESDVertex* vtxESD = 0x0;
-     if      (fVtxTPC) vtxESD = fESD->GetPrimaryVertexTPC() ;
-     else if (fVtxSPD) vtxESD = fESD->GetPrimaryVertexSPD() ;
-     else              vtxESD = fESD->GetPrimaryVertexTracks() ;
-     if(!vtxESD){
-       return;
-     }
-     Int_t nCtrb = vtxESD->GetNContributors();
-     if(nCtrb < fMinNbContributors) return;
-     Double_t zPosition = vtxESD->GetZ();
-     if(TMath::Abs(zPosition) > fRangePrimaryVertexZ) return;     
-
-   }
-
-   //printf("Primary vertex passed\n");
-
-   fNEvents->Fill(1);
-
-   // In total
-   Int_t nbTrdTracks = 0;
-   // standalone
-   Int_t nbTrdTracksStandalone = 0;
-   // offline
-   Int_t nbTrdTracksOffline = 0;
-   // TPC
-   Int_t nbtrackTPC = 0;
-
-   Double_t nbTracks = fESD->GetNumberOfTracks();
-   //printf("Number of tracks %f\n",nbTracks);  
-
-   if (nbTracks <= 0.0) {
-     
-     if(fDebug > 1) {
-       fNbTRDTrack->Fill(nbTrdTracks);
-       fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);
-       fNbTRDTrackOffline->Fill(nbTrdTracksOffline);
-     }
-     PostData(0, fListHist);
-     return;
-   }
-
-   //fESDfriend = (AliESDfriend *)fESD->FindListObject("AliESDfriend");
-   fESDfriend = (AliESDfriend *)fESD->FindListObject("AliESDfriend");
-   fESD->SetESDfriend(fESDfriend);
-   //if(!fESDfriend) return;   
-   
-   //printf("has friends\n");
+  //  AliLog::SetGlobalLogLevel(AliLog::kError);
+  //  cout << "AliTRDCalibTask::Exec() 1" << endl;
+  fESD = dynamic_cast<AliESDEvent*>(fInputEvent);
+  if(!fESD){
+    AliError("ESD Event missing");
+    PostData(1, fListHist);
+    return;
+  }
+  
+  //printf("Counter %d\n",fCounter);
+  
+  fCounter++;
+  //cout << "maxEvent = " << fMaxEvent << endl;
+  //if(fCounter%100==0) cout << "fCounter = " << fCounter << endl;
+  if((fMaxEvent != 0) && (fMaxEvent < fCounter)) return;
+  //if(fCounter%100==0) cout << "fCounter1 = " << fCounter << endl;
+  //cout << "event = " << fCounter << endl;
+  
+  //printf("Counter %d\n",fCounter);
+  
+  ///////////////////
+  // Check trigger
+  ///////////////////
+  Bool_t pass = kTRUE;
+  Int_t numberOfTriggerSelected = fSelectedTrigger->GetEntriesFast();
+  //printf("numberofTriggerSelected %d\n",numberOfTriggerSelected);
+  if(fRejected) {
+    pass = kTRUE;
+    for(Int_t k = 0; k < numberOfTriggerSelected; k++){
+      const TObjString *const obString=(TObjString*)fSelectedTrigger->At(k);
+      const TString tString=obString->GetString();
+      if(fESD->IsTriggerClassFired((const char*)tString)) {
+       pass = kFALSE;
+      }
+    }
+  }
+  else {
+    pass = kFALSE;
+    for(Int_t k = 0; k < numberOfTriggerSelected; k++){
+      const TObjString *const obString=(TObjString*)fSelectedTrigger->At(k);
+      const TString tString=obString->GetString();
+      if(fESD->IsTriggerClassFired((const char*)tString)) {
+       pass = kTRUE;
+      }
+    }
+  }
+  if(!pass) {
+    PostData(1, fListHist);
+    return;
+  }   
+  //printf("Class Fired %s\n",(const char*)fESD->GetFiredTriggerClasses());
+  //printf("Trigger passed\n");
+  
+  ///////////////////////////////
+  // Require a primary vertex
+  //////////////////////////////
+  if(fRequirePrimaryVertex) {
+    const AliESDVertex* vtxESD = 0x0;
+    if      (fVtxTPC) vtxESD = fESD->GetPrimaryVertexTPC() ;
+    else if (fVtxSPD) vtxESD = fESD->GetPrimaryVertexSPD() ;
+    else              vtxESD = fESD->GetPrimaryVertexTracks() ;
+    if(!vtxESD){
+      PostData(1, fListHist);
+      return;
+    }
+    Int_t nCtrb = vtxESD->GetNContributors();
+    if(nCtrb < fMinNbContributors) {
+      PostData(1, fListHist);     
+      return;
+    }
+    Double_t zPosition = vtxESD->GetZ();
+    if(TMath::Abs(zPosition) > fRangePrimaryVertexZ) {
+      PostData(1, fListHist);
+      return;
+    }     
+    
+  }
+  
+  //printf("Primary vertex passed\n");
+  
+  fNEvents->Fill(1);
+  
+  // In total
+  Int_t nbTrdTracks = 0;
+  // standalone
+  Int_t nbTrdTracksStandalone = 0;
+  // offline
+  Int_t nbTrdTracksOffline = 0;
+  // TPC
+  Int_t nbtrackTPC = 0;
+  
+  Double_t nbTracks = fESD->GetNumberOfTracks();
+  //printf("Number of tracks %f\n",nbTracks);  
+  
+  if (nbTracks <= 0.0) {
+    
+    if(fDebug > 1) {
+      fNbTRDTrack->Fill(nbTrdTracks);
+      fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);
+      fNbTRDTrackOffline->Fill(nbTrdTracksOffline);
+    }
+    PostData(1, fListHist);
+    return;
+  }
+  
+  fESDfriend = dynamic_cast<AliESDfriend*> (fESD->FindListObject("AliESDfriend"));
+  if(!fESDfriend){
+    AliError("fESDfriend not available");
+    PostData(1, fListHist);
+    return;
+  }
+  
+  //printf("has friends\n");
 
-   ////////////////////////////////////
+  /*
+  ////////////////////////////////////
    // Check the number of TPC tracks
    ///////////////////////////////////
    //printf("Nb of tracks %f\n",nbTracks);
@@ -519,31 +545,16 @@ void AliTRDCalibTask::CreateOutputObjects()
      fkEsdTrack = fESD->GetTrack(itrk);
      ULong_t status = fkEsdTrack->GetStatus(); 
      if(status&(AliESDtrack::kTPCout)) nbtrackTPC++;
-     // Check that the calibration object is here
-     if(fESDfriend && (fESDfriend->GetTrack(itrk))) {
-       fFriendTrack = new AliESDfriendTrack(*(fESDfriend->GetTrack(itrk)));
-       Int_t counteer = 0;
-       Int_t icalib=0;
-       //printf("Found friend track %d\n",itrk);
-       while((fCalibObject = (TObject *)(fFriendTrack->GetCalibObject(icalib++)))){
-        //printf("Name %s\n",fCalibObject->IsA()->GetName());
-                if(strcmp(fCalibObject->IsA()->GetName(), "AliTRDtrackV1") != 0) continue;
-        counteer++;
-       }
-       //printf("TRDntracklets %d, TRDntrackletsPID %d\n",fkEsdTrack->GetTRDntracklets(),fkEsdTrack->GetTRDntrackletsPID());
-       if(counteer > 0) {
-        nbTrdTracks++;      
-        if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {
-          nbTrdTracksStandalone++;
-        }
-        if((status&(AliESDtrack::kTRDin))) {
-          nbTrdTracksOffline++;
-        }
-       }
-       if(fFriendTrack) delete fFriendTrack;
+     if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {
+       nbTrdTracks++;    
+       nbTrdTracksStandalone++;
+     }
+     if((status&(AliESDtrack::kTRDin))) {
+       nbTrdTracks++;    
+       nbTrdTracksOffline++;
      }
    }
-   //printf("Number of TPC tracks %d, TRD %d\n",nbtrackTPC,nbTrdTracks);
+   
    if((nbtrackTPC>0) && (nbTrdTracks > (3.0*nbtrackTPC))) pass = kFALSE;
    
    if(fDebug > 1) {
@@ -556,200 +567,229 @@ void AliTRDCalibTask::CreateOutputObjects()
    }
 
    if(!pass) {
-     PostData(0, fListHist);
+     PostData(1, fListHist);
      return;
    }
-   //printf("Pass\n");  
-
-   /////////////////////////////////////
-   // Loop on AliESDtrack
-   ////////////////////////////////////
-   //printf("Nb of tracks %f\n",nbTracks);      
-   for(int itrk=0; itrk < nbTracks; itrk++){
-
-     // Get ESD track
-     fkEsdTrack = fESD->GetTrack(itrk);
-     if(!fkEsdTrack) continue;
-
-     // Quality cuts on the AliESDtrack
-     if((fEsdTrackCuts) && (!fEsdTrackCuts->IsSelected((AliVParticle *)fkEsdTrack))) {
-       //printf("Not a good track\n");
-       continue;
-     }
-
-     // First Absolute gain calibration
-     Int_t trdNTracklets = (Int_t) fkEsdTrack->GetTRDntracklets();
-     Int_t trdNTrackletsPID = (Int_t) fkEsdTrack->GetTRDntrackletsPID(); 
-     if((trdNTracklets > 0) && (trdNTrackletsPID > 0)) {
-       for(Int_t iPlane = 0; iPlane < 6; iPlane++){
-        Double_t slide = fkEsdTrack->GetTRDslice(iPlane);
-        //printf("Number of slide %d\n",fkEsdTrack->GetNumberOfTRDslices());
-        Double_t momentum = fkEsdTrack->GetTRDmomentum(iPlane);
-        //printf("momentum %f, slide %f\n",momentum,slide);
-        if(slide > 0.0) fAbsoluteGain->Fill(slide*8.0/100.0,momentum); 
-       }
-     }     
-     
-     // Other cuts
-     Bool_t good = kTRUE;
-     Bool_t standalonetrack = kFALSE;
-     Bool_t offlinetrack = kFALSE;
-     ULong_t status = fkEsdTrack->GetStatus();
-     
-     if(!(fESDfriend->GetTrack(itrk)))  continue;   
-     
-     fFriendTrack = new AliESDfriendTrack(*(fESDfriend->GetTrack(itrk)));
-     
-     //////////////////////////////////////
-     // Loop on calibration objects
-     //////////////////////////////////////
-     Int_t icalib=0;
-     while((fCalibObject = (TObject *)(fFriendTrack->GetCalibObject(icalib++)))){
-       if(strcmp(fCalibObject->IsA()->GetName(), "AliTRDtrackV1") != 0) continue;
-       //printf("Find the calibration object\n");
-
-       if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {
-        standalonetrack = kTRUE;
-       }
-       if((status&(AliESDtrack::kTRDin))) {
-        offlinetrack = kTRUE;
-       }
-       if(fOfflineTracks){
-        if(!offlinetrack){
-          good = kFALSE;
-        }
-       }
-       else if(fStandaloneTracks){
-        if(!standalonetrack){
-          good = kFALSE;
-        }
-       }
-       
-       fTrdTrack = (AliTRDtrackV1 *)fCalibObject;
-       if(good) {
-        fTRDCalibraFillHisto->UpdateHistogramsV1(fTrdTrack);
-        //printf("Fill fTRDCalibraFillHisto\n");
-       }
-       
-       //////////////////////////////////
-       // Debug 
-       ////////////////////////////////
-
-       if(fDebug > 0) {
-        
-        //printf("Enter debug\n");
-        
-        Int_t nbtracklets = 0;
-        
-        // Check some stuff
-        Bool_t standalonetracklet = kFALSE;  
-        const AliTRDseedV1 *tracklet = 0x0;
-        //////////////////////////////////////
-        // Loop tracklets
-        ///////////////////////////////////// 
-        for(Int_t itr = 0; itr < 6; itr++){
-          
-          if(!(tracklet = fTrdTrack->GetTracklet(itr))) continue;
-          if(!tracklet->IsOK()) continue;
-          nbtracklets++;
-          standalonetracklet = kFALSE; 
-          if(tracklet->IsStandAlone()) standalonetracklet = kTRUE;
-
-          Int_t nbclusters = 0;
-          Double_t phtb[AliTRDseedV1::kNtb];
-          memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
-          Double_t sum = 0.0;
-          Float_t normalisation = 6.67;
-          Int_t detector = 0;
-          Int_t sector = 0;
-          //Int_t crossrow = 0;
-          
-          // Check no shared clusters
-          //for(int icc=AliTRDseedV1::kNtb; icc<AliTRDseedV1::kNclusters; icc++){
-          //  if((fcl = tracklet->GetClusters(icc)))  crossrow = 1;
-          // }
-          
-          // Loop on clusters
-          for(int ic=0; ic<AliTRDseedV1::kNtb; ic++){
-            
-            if(!(fCl = tracklet->GetClusters(ic))) continue;
-            nbclusters++;
-            Int_t time = fCl->GetPadTime();
-            Float_t ch =  tracklet->GetdQdl(ic);
-            Float_t qcl = TMath::Abs(fCl->GetQ());
-            detector = fCl->GetDetector();       
-            // Add the charge if shared cluster
-            if((ic+AliTRDseedV1::kNtb) < AliTRDseedV1::kNclusters) {
-              if((fCl = tracklet->GetClusters(ic+AliTRDseedV1::kNtb))) {
-                qcl += TMath::Abs(fCl->GetQ());
-                //printf("Add the cluster charge\n");
-              }
-            }
-            if((time>-1) && (time<fNbTimeBins)) phtb[time]=qcl;
-            if((fCalDetGain) && (fCalDetGain->GetValue(detector) > 0.0)) sum += ch*fCalDetGain->GetValue(detector)/normalisation;      
-            else sum += ch/normalisation;
-            
-            if(fDebug > 1) {
-              fNbTimeBin->Fill(time);
-              if(tracklet->IsStandAlone()) fNbTimeBinStandalone->Fill(time);
-              else fNbTimeBinOffline->Fill(time);
-            }
-          }
-          sector = AliTRDgeometry::GetSector(detector);
-          
-          if(fDebug > 1) {
-            fNbTracklets->Fill(detector);
-            if(tracklet->IsStandAlone()) fNbTrackletsStandalone->Fill(detector);
-            else fNbTrackletsOffline->Fill(detector);
-          
-            fNbClusters->Fill(nbclusters);
-            if(tracklet->IsStandAlone())  fNbClustersStandalone->Fill(nbclusters);
-            else  fNbClustersOffline->Fill(nbclusters);
-          }       
-
-          if(fDebug > 0) {
-            if((nbclusters > fLow) && (nbclusters < fHigh)){
-              if(fRelativeScale > 0.0) sum = sum/fRelativeScale;              
-              fCH2dSM->Fill(sum,sector+0.5);
-              fCH2dSum->Fill(sum,0.5);
-              Bool_t checknoise = kTRUE;
-              if(fMaxCluster > 0) {
-                if(phtb[0] > fMaxCluster) checknoise = kFALSE;
-                if(fNbTimeBins > fNbMaxCluster) {
-                  for(Int_t k = (fNbTimeBins-fNbMaxCluster); k < fNbTimeBins; k++){
-                    if(phtb[k] > fMaxCluster) checknoise = kFALSE;
-                  }
-                }
-              }
-              if(checknoise) {        
-                for(int ic=0; ic<fNbTimeBins; ic++){
-                  if(fFillZero) {
-                    fPH2dSum->Fill((Double_t)(ic/10.0),0.5,(Double_t)phtb[ic]);
-                    fPH2dSM->Fill((Double_t)(ic/10.0),sector+0.5,(Double_t)phtb[ic]);
-                  }
-                  else {
-                    if(phtb[ic] > 0.0) {
-                      fPH2dSum->Fill((Double_t)(ic/10.0),0.0,(Double_t)phtb[ic]);
-                      fPH2dSM->Fill((Double_t)(ic/10.0),sector+0.5,(Double_t)phtb[ic]);
-                    }
-                  }
-                }
-              }
-            }
-          }
-        } // loop on tracklets
+  */
+  
+  /////////////////////////////////////
+  // Loop on AliESDtrack
+  ////////////////////////////////////
+  //printf("Nb of tracks %f\n",nbTracks);      
+  for(int itrk=0; itrk < nbTracks; ++itrk){
     
-       } // debug
-       
-     }// while calibration objects
-     
-     delete fFriendTrack;
+    // Get ESD track
+    fkEsdTrack = fESD->GetTrack(itrk);
+    if(!fkEsdTrack) continue;
+    ULong_t status = fkEsdTrack->GetStatus(); 
+    if(status&(AliESDtrack::kTPCout)) ++nbtrackTPC;
+    
+    // Quality cuts on the AliESDtrack
+    if((fEsdTrackCuts) && (!fEsdTrackCuts->IsSelected((AliVParticle *)fkEsdTrack))) {
+      //printf("Not a good track\n");
+      continue;
+    }
+    
+    // First Absolute gain calibration
+    Int_t trdNTracklets = (Int_t) fkEsdTrack->GetTRDntracklets();
+    Int_t trdNTrackletsPID = (Int_t) fkEsdTrack->GetTRDntrackletsPID(); 
+    if((trdNTracklets > 0) && (trdNTrackletsPID > 0)) {
+      for(Int_t iPlane = 0; iPlane < 6; ++iPlane){
+       //Double_t slide = fkEsdTrack->GetTRDslice(iPlane);
+       //printf("Number of slide %d\n",fkEsdTrack->GetNumberOfTRDslices());
+       //Double_t momentum = fkEsdTrack->GetTRDmomentum(iPlane);
+       //printf("momentum %f, slide %f\n",momentum,slide);
+       if(fkEsdTrack->GetTRDslice(iPlane) > 0.0) 
+         fAbsoluteGain->Fill(fkEsdTrack->GetTRDslice(iPlane)*8.0/100.0,
+                             fkEsdTrack->GetTRDmomentum(iPlane)); 
+      }
+    }     
+    
+    // Other cuts
+    Bool_t good = kTRUE;
+    Bool_t standalonetrack = kFALSE;
+    Bool_t offlinetrack = kFALSE;
+    //ULong_t status = fkEsdTrack->GetStatus();
+    
+    fFriendTrack = fESDfriend->GetTrack(itrk);
+    if(!fFriendTrack)  continue;
+    //////////////////////////////////////
+    // Loop on calibration objects
+    //////////////////////////////////////
+    Int_t icalib=0;
+    Int_t nTRDtrackV1=0;
+    while((fCalibObject = (TObject *)(fFriendTrack->GetCalibObject(icalib++)))){
+      if(strcmp(fCalibObject->IsA()->GetName(), "AliTRDtrackV1") != 0) continue;
+      //printf("Find the calibration object\n");
+      ++nTRDtrackV1;
+      
+      if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {
+       standalonetrack = kTRUE;
+      }
+      if((status&(AliESDtrack::kTRDin))) {
+       offlinetrack = kTRUE;
+      }
+      if(fOfflineTracks){
+       if(!offlinetrack){
+         good = kFALSE;
+       }
+      }
+      else if(fStandaloneTracks){
+       if(!standalonetrack){
+         good = kFALSE;
+       }
+      }
+      
+      fTrdTrack = (AliTRDtrackV1 *)fCalibObject;
+      if(good) {
+       //cout << "good" << endl;
+       fTRDCalibraFillHisto->UpdateHistogramsV1(fTrdTrack);
+       //printf("Fill fTRDCalibraFillHisto\n");
+      }
+      
+      //////////////////////////////////
+      // Debug 
+      ////////////////////////////////
+      
+      if(fDebug > 0) {
+       
+       //printf("Enter debug\n");
+       
+       Int_t nbtracklets = 0;
+       
+       // Check some stuff
+       Bool_t standalonetracklet = kFALSE;  
+       const AliTRDseedV1 *tracklet = 0x0;
+       //////////////////////////////////////
+       // Loop tracklets
+       ///////////////////////////////////// 
+       Int_t nbclusters=0;
+       Double_t phtb[AliTRDseedV1::kNtb];
+        memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
+       Double_t sum = 0.0;
+       Float_t normalisation = 6.67;
+       Int_t detector = 0;
+       Int_t sector = 0;
+       for(Int_t itr = 0; itr < 6; ++itr){
+         
+         if(!(tracklet = fTrdTrack->GetTracklet(itr))) continue;
+         if(!tracklet->IsOK()) continue;
+         ++nbtracklets;
+         standalonetracklet = kFALSE; 
+         if(tracklet->IsStandAlone()) standalonetracklet = kTRUE;
+         
+         nbclusters = 0;
+         memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));
+         sum = 0.0;
+         detector = 0;
+         sector = 0;
+         //Int_t crossrow = 0;
+         
+         // Check no shared clusters
+         //for(int icc=AliTRDseedV1::kNtb; icc<AliTRDseedV1::kNclusters; icc++){
+         //  if((fcl = tracklet->GetClusters(icc)))  crossrow = 1;
+         // }
+         
+         // Loop on clusters
+         Int_t time = 0;
+         Float_t ch = 0;
+         Float_t qcl = 0;
+         for(int ic=0; ic<AliTRDseedV1::kNtb; ++ic){
+           
+           if(!(fCl = tracklet->GetClusters(ic))) continue;
+           ++nbclusters;
+           time = fCl->GetPadTime();
+           ch =  tracklet->GetdQdl(ic);
+           qcl = TMath::Abs(fCl->GetQ());
+           detector = fCl->GetDetector();        
+           // Add the charge if shared cluster
+           if((ic+AliTRDseedV1::kNtb) < AliTRDseedV1::kNclusters) {
+             if((fCl = tracklet->GetClusters(ic+AliTRDseedV1::kNtb))) {
+               qcl += TMath::Abs(fCl->GetQ());
+               //printf("Add the cluster charge\n");
+             }
+           }
+           if((time>-1) && (time<fNbTimeBins)) phtb[time]=qcl;
+           if((fCalDetGain) && (fCalDetGain->GetValue(detector) > 0.0)) sum += ch*fCalDetGain->GetValue(detector)/normalisation;       
+           else sum += ch/normalisation;
+           
+           if(fDebug > 1) {
+             fNbTimeBin->Fill(time);
+             if(tracklet->IsStandAlone()) fNbTimeBinStandalone->Fill(time);
+             else fNbTimeBinOffline->Fill(time);
+           }
+         }
+         sector = AliTRDgeometry::GetSector(detector);
+         
+         if(fDebug > 1) {
+           fNbTracklets->Fill(detector);
+           if(tracklet->IsStandAlone()) fNbTrackletsStandalone->Fill(detector);
+           else fNbTrackletsOffline->Fill(detector);
+           
+           fNbClusters->Fill(nbclusters);
+           if(tracklet->IsStandAlone())  fNbClustersStandalone->Fill(nbclusters);
+           else  fNbClustersOffline->Fill(nbclusters);
+         }        
+         
+         if(fDebug > 0) {
+           if((nbclusters > fLow) && (nbclusters < fHigh)){
+             if(fRelativeScale > 0.0) sum = sum/fRelativeScale;               
+             fCH2dSM->Fill(sum,sector+0.5);
+             fCH2dSum->Fill(sum,0.5);
+             Bool_t checknoise = kTRUE;
+             if(fMaxCluster > 0) {
+               if(phtb[0] > fMaxCluster) checknoise = kFALSE;
+               if(fNbTimeBins > fNbMaxCluster) {
+                 for(Int_t k = (fNbTimeBins-fNbMaxCluster); k < fNbTimeBins; k++){
+                   if(phtb[k] > fMaxCluster) checknoise = kFALSE;
+                 }
+               }
+             }
+             if(checknoise) {         
+               for(int ic=0; ic<fNbTimeBins; ic++){
+                 if(fFillZero) {
+                   fPH2dSum->Fill((Double_t)(ic/10.0),0.5,(Double_t)phtb[ic]);
+                   fPH2dSM->Fill((Double_t)(ic/10.0),sector+0.5,(Double_t)phtb[ic]);
+                 }
+                 else {
+                   if(phtb[ic] > 0.0) {
+                     fPH2dSum->Fill((Double_t)(ic/10.0),0.0,(Double_t)phtb[ic]);
+                     fPH2dSM->Fill((Double_t)(ic/10.0),sector+0.5,(Double_t)phtb[ic]);
+                   }
+                 }
+               }
+             }
+           }
+         }
+       } // loop on tracklets
+       
+      } // debug
+      
+    }// while calibration objects
+    if(nTRDtrackV1 > 0) {
+      ++nbTrdTracks;      
+      if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {
+       ++nbTrdTracksStandalone;
+      }
+      if((status&(AliESDtrack::kTRDin))) {
+       ++nbTrdTracksOffline;
+      }
+    }
+    //delete fFriendTrack;
+  } // loop ESD track
+  
+  if(fDebug > 1) {
+    fNbTRDTrack->Fill(nbTrdTracks);
+    fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);
+    fNbTRDTrackOffline->Fill(nbTrdTracksOffline);
+    fNbTPCTRDtrack->Fill(nbTrdTracks,nbtrackTPC);
+  }
+  
+  // Post output data
+  PostData(1, fListHist);
+  //cout << "AliTRDCalibTask::Exec() OUT" << endl;
+}
      
-   } // loop ESD track
-   
-   // Post output data
-   PostData(0, fListHist);
- }     
 //________________________________________________________________________
 void AliTRDCalibTask::Terminate(Option_t *) 
 {
@@ -1288,3 +1328,4 @@ Long64_t AliTRDCalibTask::Merge(TCollection *li) {
   return 0;
   
 }
+
index d2bf604..39a9f42 100644 (file)
@@ -28,17 +28,17 @@ class AliESDtrackCuts;
 class AliTRDCalDet;
 
 #include "TObjString.h"
-#include "AliAnalysisTask.h" 
+#include "AliAnalysisTaskSE.h" 
 #include "TMath.h"
 
-class AliTRDCalibTask : public AliAnalysisTask {
+class AliTRDCalibTask : public AliAnalysisTaskSE {
  public:
   AliTRDCalibTask(const char *name = "AliTRDCalibTask");
   virtual ~AliTRDCalibTask();
   
-  virtual void   ConnectInputData(Option_t *);
-  virtual void   CreateOutputObjects();
-  virtual void   Exec(Option_t *);
+  //  virtual void   ConnectInputData(Option_t *);
+  virtual void   UserCreateOutputObjects();
+  virtual void   UserExec(Option_t *);
   virtual void   Terminate(Option_t *);
   virtual Bool_t Load(const Char_t *filename);
   virtual Bool_t Load(TList *lister);
@@ -157,3 +157,4 @@ class AliTRDCalibTask : public AliAnalysisTask {
 };
 
 #endif
+
index 201e04e..57a240a 100644 (file)
@@ -126,6 +126,7 @@ void AliTRDCalibraFillHisto::Terminate()
 AliTRDCalibraFillHisto::AliTRDCalibraFillHisto()
   :TObject()
   ,fGeo(0)
+  ,fCalibDB(0)
   ,fIsHLT(kFALSE)
   ,fCH2dOn(kFALSE)
   ,fPH2dOn(kFALSE)
@@ -191,13 +192,14 @@ AliTRDCalibraFillHisto::AliTRDCalibraFillHisto()
   fNumberUsedPh[1]       = 0;
   
   fGeo = new AliTRDgeometry();
-
+  fCalibDB = AliTRDcalibDB::Instance();
 }
 
 //______________________________________________________________________________________
 AliTRDCalibraFillHisto::AliTRDCalibraFillHisto(const AliTRDCalibraFillHisto &c)
   :TObject(c)
   ,fGeo(0)
+  ,fCalibDB(0)
   ,fIsHLT(c.fIsHLT)
   ,fCH2dOn(c.fCH2dOn)
   ,fPH2dOn(c.fPH2dOn)
@@ -277,6 +279,7 @@ AliTRDCalibraFillHisto::AliTRDCalibraFillHisto(const AliTRDCalibraFillHisto &c)
     delete fGeo;
   }
   fGeo = new AliTRDgeometry();
+  fCalibDB = AliTRDcalibDB::Instance();
 }
 
 //____________________________________________________________________________________
@@ -649,11 +652,18 @@ Bool_t AliTRDCalibraFillHisto::UpdateHistogramsV1(const AliTRDtrackV1 *t)
   Bool_t         newtr   = kTRUE;              // new track
   
   // Get cal
-  AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
+  //  AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
+  /*
   if (!cal) {
     AliInfo("Could not get calibDB");
     return kFALSE;
   }
+*/
+  if (!fCalibDB) {
+    AliInfo("Could not get calibDB");
+    return kFALSE;
+  }
+
   
   ///////////////////////////
   // loop over the tracklet
@@ -687,9 +697,9 @@ Bool_t AliTRDCalibraFillHisto::UpdateHistogramsV1(const AliTRDtrackV1 *t)
       if( fCalROCGain ){ 
        if(!fIsHLT){    
          fCalROCGain->~AliTRDCalROC();
-         new(fCalROCGain) AliTRDCalROC(*(cal->GetGainFactorROC(detector)));
+         new(fCalROCGain) AliTRDCalROC(*(fCalibDB->GetGainFactorROC(detector)));
        }
-      }else fCalROCGain = new AliTRDCalROC(*(cal->GetGainFactorROC(detector)));
+      }else fCalROCGain = new AliTRDCalROC(*(fCalibDB->GetGainFactorROC(detector)));
       
       // reset
       fDetectorPreviousTrack = detector;
index 04354e6..1d1cd27 100644 (file)
@@ -42,6 +42,7 @@ class AliTRDseedV1;
 class AliTRDgeometry;
 class AliTRDCalDet;
 class AliTRDCalROC;
+class AliTRDcalibDB;
 
 class AliTRDrawFastStream;
 class AliTRDdigitsManager;
@@ -173,6 +174,8 @@ AliTRDCalibraVector *GetCalibraVector() const                                { r
 
   // Geometry
   AliTRDgeometry  *fGeo;                    //! The TRD geometry
+  // calibration DB
+  AliTRDcalibDB   *fCalibDB;                //! The pointer to the TRDcalibDB instance
 
   // Is HLT
           Bool_t   fIsHLT;                  // Now if HLT, the per detector
index 91c622e..2cd0c3f 100644 (file)
 
 /*
   Responsible: marian.ivanov@cern.ch 
-  Code to analyze the TRD calibration and to produce OCDB entries  
+  Code to analyze the TPC calibration and to produce OCDB entries  
 
 
    .x ~/rootlogon.C
    gSystem->Load("libANALYSIS");
-   gSystem->Load("libTRDcalib");
+   gSystem->Load("libTPCcalib");
 
    AliTRDPreprocessorOffline proces;
    TString ocdbPath="local:////"
 ClassImp(AliTRDPreprocessorOffline)
 
 AliTRDPreprocessorOffline::AliTRDPreprocessorOffline():
-  TNamed("TRDPreprocessorOffline","TRDPreprocessorOffline"),
+  TNamed("TPCPreprocessorOffline","TPCPreprocessorOffline"),
   fMethodSecond(kTRUE),
+  fNameList("TRDCalib"),
+  fCalDetGainUsed(0x0),
   fCH2d(0x0),
   fPH2d(0x0),
   fPRF2d(0x0),
@@ -82,6 +84,7 @@ AliTRDPreprocessorOffline::~AliTRDPreprocessorOffline() {
   // Destructor
   //
 
+  if(fCalDetGainUsed) delete fCalDetGainUsed;
   if(fCH2d) delete fCH2d;
   if(fPH2d) delete fPH2d;
   if(fPRF2d) delete fPRF2d;
@@ -142,6 +145,7 @@ void AliTRDPreprocessorOffline::CalibGain(const Char_t* file, Int_t startRunNumb
   // 2. extraction of the information
   //
   AnalyzeGain();
+  if(fCalDetGainUsed) CorrectFromDetGainUsed();
   //
   // 3. Append QA plots
   //
@@ -192,7 +196,7 @@ Bool_t AliTRDPreprocessorOffline::ReadGainGlobal(const Char_t* fileName){
   // read calibration entries from file
   // 
   TFile fcalib(fileName);
-  TObjArray * array = (TObjArray*)fcalib.Get("TRDCalib");
+  TObjArray * array = (TObjArray*)fcalib.Get(fNameList);
   if (array){
     TH2I *ch2d = (TH2I *) array->FindObject("CH2d");
     if(!ch2d) return kFALSE;
@@ -218,7 +222,7 @@ Bool_t AliTRDPreprocessorOffline::ReadVdriftT0Global(const Char_t* fileName){
   // read calibration entries from file
   // 
   TFile fcalib(fileName);
-  TObjArray * array = (TObjArray*)fcalib.Get("TRDCalib");
+  TObjArray * array = (TObjArray*)fcalib.Get(fNameList);
   if (array){
     TProfile2D *ph2d = (TProfile2D *) array->FindObject("PH2d");
     if(!ph2d) return kFALSE;
@@ -242,7 +246,7 @@ Bool_t AliTRDPreprocessorOffline::ReadVdriftLinearFitGlobal(const Char_t* fileNa
   // read calibration entries from file
   // 
   TFile fcalib(fileName);
-  TObjArray * array = (TObjArray*)fcalib.Get("TRDCalib");
+  TObjArray * array = (TObjArray*)fcalib.Get(fNameList);
   if (array){
     fAliTRDCalibraVdriftLinearFit = (AliTRDCalibraVdriftLinearFit *) array->FindObject("AliTRDCalibraVdriftLinearFit");
     //fNEvents = (TH1I *) array->FindObject("NEvents");
@@ -263,7 +267,7 @@ Bool_t AliTRDPreprocessorOffline::ReadPRFGlobal(const Char_t* fileName){
   // read calibration entries from file
   // 
   TFile fcalib(fileName);
-  TObjArray * array = (TObjArray*)fcalib.Get("TRDCalib");
+  TObjArray * array = (TObjArray*)fcalib.Get(fNameList);
   if (array){
     TProfile2D *prf2d = (TProfile2D *) array->FindObject("PRF2d");
     if(!prf2d) return kFALSE;
@@ -343,10 +347,11 @@ Bool_t AliTRDPreprocessorOffline::AnalyzeVdriftT0(){
       (nbfit        >= 0.5*nbE) && (nbE > 30) && (nbfitSuccess > 30)) {
     //printf("Pass the cut for VdriftT0\n");
     // create the cal objects
-    calibra->RemoveOutliers(1,kTRUE);
-    calibra->PutMeanValueOtherVectorFit(1,kTRUE);
-    calibra->RemoveOutliers2(kTRUE);
-    calibra->PutMeanValueOtherVectorFit2(1,kTRUE);
+    calibra->RemoveOutliers(1,kFALSE);
+    calibra->PutMeanValueOtherVectorFit(1,kFALSE);
+    calibra->RemoveOutliers2(kFALSE);
+    calibra->PutMeanValueOtherVectorFit2(1,kFALSE);
+    //
     TObjArray object = calibra->GetVectorFit();
     AliTRDCalDet *calDetVdrift = calibra->CreateDetObjectVdrift(&object);
     TH1F *coefVdriftPH  = calDetVdrift->MakeHisto1DAsFunctionOfDet();
@@ -382,6 +387,7 @@ Bool_t AliTRDPreprocessorOffline::AnalyzeVdriftLinearFit(){
 
   AliTRDCalibraFit *calibra = AliTRDCalibraFit::Instance();
   calibra->SetMinEntries(800); // If there is less than 1000 entries in the histo: no fit
+  fAliTRDCalibraVdriftLinearFit->FillPEArray();
   calibra->AnalyseLinearFitters(fAliTRDCalibraVdriftLinearFit);
 
   //Int_t nbtg        = 540;
@@ -393,6 +399,11 @@ Bool_t AliTRDPreprocessorOffline::AnalyzeVdriftLinearFit(){
   // enough statistics
   if ((nbfit        >= 0.5*nbE) && (nbE > 30)) {
     // create the cal objects
+    //calibra->RemoveOutliers(1,kTRUE);
+    calibra->PutMeanValueOtherVectorFit(1,kTRUE);
+    //calibra->RemoveOutliers2(kTRUE);
+    calibra->PutMeanValueOtherVectorFit2(1,kTRUE);
+    //
     TObjArray object  = calibra->GetVectorFit();
     AliTRDCalDet *calDetVdrift = calibra->CreateDetObjectVdrift(&object,kTRUE);
     TH1F *coefDriftLinear      = calDetVdrift->MakeHisto1DAsFunctionOfDet();
@@ -451,6 +462,22 @@ Bool_t AliTRDPreprocessorOffline::AnalyzePRF(){
   
 }
 
+void AliTRDPreprocessorOffline::CorrectFromDetGainUsed() {
+
+  AliTRDCalDet *calDetGain = (AliTRDCalDet *) fCalibObjects->At(kGain);
+  if(!calDetGain) return;
+
+  for(Int_t det = 0; det < 540; det++) {
+    
+    Float_t gaininit = fCalDetGainUsed->GetValue(det);
+    Float_t gainout = calDetGain->GetValue(det);
+    
+    calDetGain->SetValue(det,gaininit*gainout);
+
+  }
+
+
+}
 
 void AliTRDPreprocessorOffline::UpdateOCDBGain(Int_t startRunNumber, Int_t endRunNumber, const Char_t *storagePath){
   //
@@ -573,3 +600,4 @@ void AliTRDPreprocessorOffline::UpdateOCDBPRF(Int_t startRunNumber, Int_t endRun
  
 
 }
+
index 6f796e7..e15b66e 100644 (file)
 
 #include "TNamed.h"
 class TObjArray;
+class AliTRDCalDet;
 class TH2I;
 class TProfile2D;
 class AliTRDCalibraVdriftLinearFit;
 class TH1I;
 class TH2F;
+class TString;
 
 
 class AliTRDPreprocessorOffline:public TNamed { 
@@ -37,6 +39,10 @@ public:
 
   void SetLinearFitForVdrift(Bool_t methodsecond) { fMethodSecond = methodsecond;};
   Bool_t GetLinearFitForVdrift() const { return fMethodSecond;};
+  void SetNameList(TString nameList) { fNameList = nameList;};
+  TString GetNameList() const { return fNameList;}; 
+  void SetCalDetGain(AliTRDCalDet *calDetGainUsed) {fCalDetGainUsed = calDetGainUsed;};
+  AliTRDCalDet *GetCalDetGain() const { return fCalDetGainUsed;};
 
   void CalibVdriftT0(const Char_t* file, Int_t startRunNumber, Int_t endRunNumber, TString ocdbStorage="");
   void CalibGain(const Char_t* file, Int_t startRunNumber, Int_t endRunNumber,  TString  ocdbStorage="");
@@ -52,6 +58,8 @@ public:
   Bool_t AnalyzeVdriftLinearFit(); 
   Bool_t AnalyzePRF(); 
   
+  void CorrectFromDetGainUsed();
+
   void UpdateOCDBT0(Int_t startRunNumber, Int_t endRunNumber, const char* storagePath);
   void UpdateOCDBVdrift(Int_t startRunNumber, Int_t endRunNumber, const char* storagePath);
   void UpdateOCDBGain(Int_t  startRunNumber, Int_t endRunNumber, const char* storagePath);
@@ -60,6 +68,8 @@ public:
   
 private:
   Bool_t fMethodSecond;                   // Second Method for drift velocity   
+  TString fNameList;                      // Name of the list
+  AliTRDCalDet *fCalDetGainUsed;          // CalDet used and to be corrected for
   TH2I *fCH2d;                            // Gain
   TProfile2D *fPH2d;                      // Drift velocity first method
   TProfile2D *fPRF2d;                     // PRF
@@ -79,3 +89,4 @@ private:
 };
 
 #endif
+
index caf11fa..67333d5 100644 (file)
@@ -29,7 +29,7 @@ AliAnalysisTask  *AddTaskTRDCalib(Int_t runNumber)
   /////////////////////////
   AliTRDCalibTask *calibTask = new AliTRDCalibTask();
   calibTask->SetHisto2d(kTRUE);
-  calibTask->SetVector2d(kTRUE);
+  calibTask->SetVector2d(kFALSE);
   calibTask->SetVdriftLinear(kTRUE);
   calibTask->SetNz(0,0);
   calibTask->SetNrphi(0,0);
@@ -39,7 +39,9 @@ AliAnalysisTask  *AddTaskTRDCalib(Int_t runNumber)
   calibTask->SetNrphi(0,2);
   calibTask->SetLow(0);
   calibTask->SetHigh(30);
-  calibTask->SetFillZero(kFALSE);
+  calibTask->SetFillZero(kTRUE);
+  calibTask->AddSelectedTriggerClass("CINT1B-ABCE-NOPF-ALL");
+  calibTask->SetReject(kFALSE);
   //calibTask->SetDebug(3);
   calibTask->SetNbTimeBins(30);
   //calibTask->SetMaxEvent(20);
@@ -78,7 +80,7 @@ AliAnalysisTask  *AddTaskTRDCalib(Int_t runNumber)
 
 
   mgr->ConnectInput(calibTask,0,cinput);
-  mgr->ConnectOutput(calibTask,0,coutput);
+  mgr->ConnectOutput(calibTask,1,coutput);
   return calibTask;
 
 }
@@ -97,4 +99,5 @@ const AliTRDCalDet *GetCalDet(Int_t runNumber){
 
   return calDet;
 
-}
\ No newline at end of file
+}
+