Upgrade following the 900 GeV data. A lot of extra things added - histograms etc...
authorhdalsgaa <hdalsgaa@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 12 Feb 2010 11:13:18 +0000 (11:13 +0000)
committerhdalsgaa <hdalsgaa@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 12 Feb 2010 11:13:18 +0000 (11:13 +0000)
PWG2/FORWARD/analysis/AliFMDAnaParameters.cxx
PWG2/FORWARD/analysis/AliFMDAnaParameters.h
PWG2/FORWARD/analysis/AliFMDAnalysisTaskCollector.cxx
PWG2/FORWARD/analysis/AliFMDAnalysisTaskCollector.h
PWG2/FORWARD/analysis/AliFMDAnalysisTaskDensity.cxx
PWG2/FORWARD/analysis/AliFMDAnalysisTaskDensity.h
PWG2/FORWARD/analysis/AliFMDAnalysisTaskSE.cxx
PWG2/FORWARD/analysis/AliFMDAnalysisTaskSE.h
PWG2/FORWARD/analysis/AliFMDAnalysisTaskSharing.cxx
PWG2/FORWARD/analysis/AliFMDAnalysisTaskSharing.h

index 959add7..4301875 100644 (file)
@@ -34,7 +34,8 @@
 #include <TMath.h>
 #include "AliTriggerAnalysis.h"
 #include "AliPhysicsSelection.h"
-//#include "AliBackgroundSelection.h"
+#include "AliBackgroundSelection.h"
+#include "AliMultiplicity.h"
 #include "AliESDEvent.h"
 #include "AliESDVertex.h"
 
@@ -88,19 +89,17 @@ AliFMDAnaParameters::AliFMDAnaParameters() :
   fSPDhighLimit(999999999),
   fCentralSelection(kFALSE)
 {
-  // Default constructor 
   fPhysicsSelection = new AliPhysicsSelection;
-  //AliBackgroundSelection* backgroundSelection = new AliBackgroundSelection("bg","bg");
+  AliBackgroundSelection* backgroundSelection = new AliBackgroundSelection("bg","bg");
   
-  //fPhysicsSelection->AddBackgroundIdentification(backgroundSelection);
+  fPhysicsSelection->AddBackgroundIdentification(backgroundSelection);
   //fPhysicsSelection->Initialize(104792);
   // Do not use this - it is only for IO 
   fgInstance = this;
+  // Default constructor 
 }
 //____________________________________________________________________
 char* AliFMDAnaParameters::GetPath(const char* species) {
-  // paths of correction and calibration objects
   
   char* path = "";
   
@@ -167,7 +166,6 @@ void AliFMDAnaParameters::Init(Bool_t forceReInit, UInt_t what)
 
 void AliFMDAnaParameters::InitBackground() {
   
-  // init background object
   //AliCDBEntry*   background = GetEntry(fgkBackgroundCorrection);
   
   TFile* fin = TFile::Open(GetPath(fgkBackgroundID));
@@ -183,7 +181,6 @@ void AliFMDAnaParameters::InitBackground() {
 
 void AliFMDAnaParameters::InitEnergyDists() {
   
-  //init energy dist object
   TFile* fin = TFile::Open(GetPath(fgkEnergyDistributionID));
   //AliCDBEntry*   edist = GetEntry(fgkEnergyDists);
   if (!fin) return;
@@ -198,7 +195,6 @@ void AliFMDAnaParameters::InitEnergyDists() {
 
 void AliFMDAnaParameters::InitEventSelectionEff() {
   
-  //init event selection object
   //AliCDBEntry*   background = GetEntry(fgkBackgroundCorrection);
   TFile* fin = TFile::Open(GetPath(fgkEventSelectionEffID));
                            
@@ -212,7 +208,7 @@ void AliFMDAnaParameters::InitEventSelectionEff() {
 
 void AliFMDAnaParameters::PrintStatus() const
 {
-  // print status of object
+  
   TString energystring;
   switch(fEnergy) {
   case k900:
@@ -283,7 +279,6 @@ void AliFMDAnaParameters::InitSharingEff() {
 }
 //____________________________________________________________________
 Float_t AliFMDAnaParameters::GetVtxCutZ() {
-  //get z vtx cut
   
   if(!fIsInit) {
     AliWarning("Not initialized yet. Call Init() to remedy");
@@ -295,7 +290,7 @@ Float_t AliFMDAnaParameters::GetVtxCutZ() {
 
 //____________________________________________________________________
 Int_t AliFMDAnaParameters::GetNvtxBins() {
-  // get number of vtx bins in analysis
+  
   if(!fIsInit) {
     AliWarning("Not initialized yet. Call Init() to remedy");
     return -1;
@@ -305,22 +300,22 @@ Int_t AliFMDAnaParameters::GetNvtxBins() {
 }
 //____________________________________________________________________
 TH1F* AliFMDAnaParameters::GetEnergyDistribution(Int_t det, Char_t ring, Float_t eta) {
-  // get energy distribution
+  
   return fEnergyDistribution->GetEnergyDistribution(det, ring, eta);
 }
 //____________________________________________________________________
 TH1F* AliFMDAnaParameters::GetEmptyEnergyDistribution(Int_t det, Char_t ring) {
-  //get energy distribution for empty events
+  
   return fEnergyDistribution->GetEmptyEnergyDistribution(det, ring);
 }
 //____________________________________________________________________
 TH1F* AliFMDAnaParameters::GetRingEnergyDistribution(Int_t det, Char_t ring) {
-  //get the energy dist for a ring
+  
   return fEnergyDistribution->GetRingEnergyDistribution(det, ring);
 }
 //____________________________________________________________________
 Float_t AliFMDAnaParameters::GetSigma(Int_t det, Char_t ring, Float_t eta) {
-  //get sigma of energy dist
+  
   if(!fIsInit) {
     AliWarning("Not initialized yet. Call Init() to remedy");
     return 0;
@@ -339,7 +334,7 @@ Float_t AliFMDAnaParameters::GetSigma(Int_t det, Char_t ring, Float_t eta) {
 
 //____________________________________________________________________
 Float_t AliFMDAnaParameters::GetMPV(Int_t det, Char_t ring, Float_t eta) {
-  //Get MPV of energy dist
+  
   if(!fIsInit) {
     AliWarning("Not initialized yet. Call Init() to remedy");
     return 0;
@@ -358,7 +353,7 @@ Float_t AliFMDAnaParameters::GetMPV(Int_t det, Char_t ring, Float_t eta) {
 }
 //____________________________________________________________________
 Float_t AliFMDAnaParameters::GetConstant(Int_t det, Char_t ring, Float_t eta) {
-  //Get constant of energy dist
+  
   if(!fIsInit) {
     AliWarning("Not initialized yet. Call Init() to remedy");
     return 0;
@@ -376,7 +371,7 @@ Float_t AliFMDAnaParameters::GetConstant(Int_t det, Char_t ring, Float_t eta) {
 }
 //____________________________________________________________________
 Float_t AliFMDAnaParameters::Get2MIPWeight(Int_t det, Char_t ring, Float_t eta) {
-  //Get 2 mip weight of energy dist
+  
   if(!fIsInit) {
     AliWarning("Not initialized yet. Call Init() to remedy");
     return 0;
@@ -396,7 +391,7 @@ Float_t AliFMDAnaParameters::Get2MIPWeight(Int_t det, Char_t ring, Float_t eta)
 }
 //____________________________________________________________________
 Float_t AliFMDAnaParameters::Get3MIPWeight(Int_t det, Char_t ring, Float_t eta) {
-  //Get 3 mip weight of energy dist
+  
   if(!fIsInit) {
     AliWarning("Not initialized yet. Call Init() to remedy");
     return 0;
@@ -446,7 +441,7 @@ Int_t AliFMDAnaParameters::GetEtaBin(Float_t eta) {
 TH2F* AliFMDAnaParameters::GetBackgroundCorrection(Int_t det, 
                                                   Char_t ring, 
                                                   Int_t vtxbin) {
-  //get background correction
+  
   if(!fIsInit) {
     AliWarning("Not initialized yet. Call Init() to remedy");
     return 0;
@@ -465,7 +460,7 @@ TH2F* AliFMDAnaParameters::GetBackgroundCorrection(Int_t det,
 
 TH1F* AliFMDAnaParameters::GetDoubleHitCorrection(Int_t det, 
                                                  Char_t ring) {
-  // get correction for double hits in p+p
+  
   if(!fIsInit) {
     AliWarning("Not initialized yet. Call Init() to remedy");
     return 0;
@@ -647,7 +642,7 @@ Bool_t AliFMDAnaParameters::IsEventTriggered(const AliESDEvent *esd) const {
   ULong64_t triggerMask = esd->GetTriggerMask();
   
   TString triggers = esd->GetFiredTriggerClasses();
-  // std::cout<<triggers.Data()<<std::endl;
+  
   //if(triggers.Contains("CINT1B-ABCE-NOPF-ALL") || triggers.Contains("CINT1B-E-NOPF-ALL")) return kTRUE;
   //else return kFALSE;
   //if(triggers.Contains("CINT1B-E-NOPF-ALL"))    return kFALSE;
@@ -667,19 +662,15 @@ Bool_t AliFMDAnaParameters::IsEventTriggered(const AliESDEvent *esd) const {
   
   switch (fTrigger) {
   case kMB1: {
-    //if (triggerMask & spdFO || ((triggerMask & v0left) || (triggerMask & v0right)))
-    //  if(/*tAna.IsOfflineTriggerFired(esd,AliTriggerAnalysis::kMB1) &&*/ (triggers.Contains("CINT1B-ABCE-NOPF-ALL")))// || triggers.Contains("CINT1-E-NOPF-ALL")))
     if(fRealData) {
       if( fPhysicsSelection->IsCollisionCandidate(esd))
        return kTRUE;
-       
-       //  triggers.Contains("CINT1C-ABCE-NOPF-ALL") || triggers.Contains("CINT1A-ABCE-NOPF-ALL"))
-    
     }
-    else
-      if(triggerMask & spdFO || ((triggerMask & v0left) || (triggerMask & v0right)))
+    else {
+      if (triggerMask & spdFO || ((triggerMask & v0left) || (triggerMask & v0right)))
        return kTRUE;
-    break;
+      break;
+    }
   }
   case kMB2: { 
     if (triggerMask & spdFO && ((triggerMask & v0left) || (triggerMask & v0right)))
@@ -696,15 +687,35 @@ Bool_t AliFMDAnaParameters::IsEventTriggered(const AliESDEvent *esd) const {
     break;
   }
   case kEMPTY: {
-    if(triggers.Contains("CBEAMB-ABCE-NOPF-ALL"))
+     const AliMultiplicity* testmult = esd->GetMultiplicity();
+    Int_t nTrackLets = testmult->GetNumberOfTracklets();
+    Int_t nClusters  = testmult->GetNumberOfSingleClusters();
+    Int_t fastOR = tAna.SPDFiredChips(esd, 0);
+    
+    const AliESDVertex* vertex = 0;
+    vertex = esd->GetPrimaryVertexSPD();
+    Bool_t vtxstatus = vertex->GetStatus();
+    if(vertex->GetNContributors() <= 0)
+      vtxstatus = kFALSE;
+    
+    if(vertex->GetZRes() > 0.1 ) 
+      vtxstatus = kFALSE;
+    
+    Bool_t v0ABG = tAna.IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0ABG);
+    
+    Bool_t v0CBG = tAna.IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0CBG);
+    Bool_t v0A = tAna.IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0A);
+    Bool_t v0C = tAna.IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0C);
+    if(triggers.Contains("CINT1A-ABCE-NOPF-ALL") || triggers.Contains("CINT1C-ABCE-NOPF-ALL")) 
       return kTRUE;
     break;
   }
+    
   }//switch
   
   return kFALSE;
 }
-
+    
 //____________________________________________________________________
 Float_t 
 AliFMDAnaParameters::GetStripLength(Char_t ring, UShort_t strip)  
index 0a1cbdb..db52e63 100644 (file)
@@ -125,15 +125,15 @@ public:
   void     SetCollisionSystem(Species collsystem) {fSpecies = collsystem;}
   void     PrintStatus() const;
   void     Print(Option_t* /* option */) const { PrintStatus(); }
-  char*    GetDndetaAnalysisName() const {return "PWG2forwardDnDeta";}
+  char*    GetDndetaAnalysisName() {return "PWG2forwardDnDeta";}
   TH1F*    GetEnergyDistribution(Int_t det, Char_t ring, Float_t eta);
   TH1F*    GetEmptyEnergyDistribution(Int_t det, Char_t ring);
   TH1F*    GetRingEnergyDistribution(Int_t det, Char_t ring);
-  AliPhysicsSelection* GetPhysicsSelection() const {return fPhysicsSelection;}
-  Bool_t   IsRealData() const {return fRealData; }
+  AliPhysicsSelection* GetPhysicsSelection() {return fPhysicsSelection;}
+  Bool_t   IsRealData() {return fRealData; }
   void     SetRealData(Bool_t realdata) {fRealData = realdata;}
-  Float_t  GetLowSPDLimit() const {return fSPDlowLimit;}
-  Float_t  GetHighSPDLimit() const {return fSPDhighLimit;}
+  Float_t  GetLowSPDLimit() {return fSPDlowLimit;}
+  Float_t  GetHighSPDLimit() {return fSPDhighLimit;}
   void     SetLowSPDLimit(Float_t cut) {fSPDlowLimit = cut;}
   void     SetHighSPDLimit(Float_t cut) {fSPDhighLimit = cut;}
   void     SetCentralTriggerSelection(Bool_t selection) {fCentralSelection = selection;}
@@ -184,16 +184,16 @@ protected:
   Bool_t fIsInit;                      //Have we been init ?
   //TObjArray*  fBackgroundArray;
   // TObjArray*  fEdistArray;
-  AliFMDAnaCalibBackgroundCorrection*         fBackground;   // background map object
-  AliFMDAnaCalibEnergyDistribution*           fEnergyDistribution; // energy distribution object
-  AliFMDAnaCalibEventSelectionEfficiency*     fEventSelectionEfficiency; // event selection distribution object
-  AliFMDAnaCalibSharingEfficiency*            fSharingEfficiency;  // sharing distribution object
+  AliFMDAnaCalibBackgroundCorrection*         fBackground;   
+  AliFMDAnaCalibEnergyDistribution*           fEnergyDistribution;
+  AliFMDAnaCalibEventSelectionEfficiency*     fEventSelectionEfficiency;
+  AliFMDAnaCalibSharingEfficiency*            fSharingEfficiency;
   //static const char* fgkBackgroundCorrection;
   //static const char* fgkEnergyDists;
-  static const char* fgkBackgroundID;  // background object
-  static const char* fgkEnergyDistributionID ;  //energy dist object
-  static const char* fgkEventSelectionEffID ;   // event selection object
-  static const char* fgkSharingEffID ;          // sharing efficiency object
+  static const char* fgkBackgroundID;
+  static const char* fgkEnergyDistributionID ;
+  static const char* fgkEventSelectionEffID ;
+  static const char* fgkSharingEffID ;
   
   TVector2 fCorner1;                  //First corner of hybrid
   TVector2 fCorner2;                  //Second corner of hybrid
@@ -207,7 +207,7 @@ protected:
   Energy   fEnergy;                   // CM energy
   MagField fMagField;                 //Magnetic field
   Species  fSpecies;                  //PbPb or pp ?
-  AliPhysicsSelection* fPhysicsSelection; // the physics selection
+  AliPhysicsSelection* fPhysicsSelection; 
   Bool_t   fRealData;                 // real or simulated
   Float_t  fSPDlowLimit ;             // low limit of SPD tracklets
   Float_t  fSPDhighLimit ;             // high limit of SPD tracklets
index 832e1ed..eeba817 100644 (file)
@@ -68,7 +68,9 @@ AliFMDAnalysisTaskCollector::AliFMDAnalysisTaskCollector()
   fOutputList(0),
   fArray(0),
   fZvtxDist(0),
-  fEvents(0)
+  fEvents(0),
+  fEmptyEventsAside(0),
+  fEmptyEventsCside(0)
 {
   // Default constructor
   
@@ -81,7 +83,9 @@ AliFMDAnalysisTaskCollector::AliFMDAnalysisTaskCollector(const char* name):
     fOutputList(0),
     fArray(0),
     fZvtxDist(0),
-    fEvents(0)
+    fEvents(0),
+    fEmptyEventsAside(0),
+    fEmptyEventsCside(0)
 {
   // Default constructor
   
@@ -159,7 +163,10 @@ void AliFMDAnalysisTaskCollector::UserExec(Option_t */*option*/)
   //Bool_t trigger = pars->IsEventTriggered(esd);
   
   Bool_t physics = pars->IsEventTriggered(esd);
-  Bool_t empty   = pars->IsEventTriggered(esd,AliFMDAnaParameters::kEMPTY);
+  //Bool_t empty   = pars->IsEventTriggered(esd,AliFMDAnaParameters::kEMPTY);
+  Bool_t emptyAside = triggers.Contains("CINT1A-ABCE-NOPF-ALL");
+  Bool_t emptyCside = triggers.Contains("CINT1C-ABCE-NOPF-ALL");
+  
   //std::cout<<physics<<"   "<<empty<<std::endl;
   //if(!trigger)
   //  physics = kFALSE;
@@ -179,10 +186,12 @@ void AliFMDAnalysisTaskCollector::UserExec(Option_t */*option*/)
   if (!fmd) return;
   if(physics)
     fEvents++;
-  else if(empty)
-    fEmptyEvents++;
+  else if(emptyAside) 
+    fEmptyEventsAside++;
+  else if(emptyCside) 
+    fEmptyEventsCside++;
   
-  if(!physics && !empty)
+  if(!physics && !emptyAside && !emptyCside)
     return;
   TH1F* Edist = 0;
   TH1F* emptyEdist = 0;
@@ -213,14 +222,27 @@ void AliFMDAnalysisTaskCollector::UserExec(Option_t */*option*/)
            Edist->Fill(mult);
            ringEdist->Fill(mult);
          }
-         else if(empty) {
+         if(emptyAside && det == 3 /*&& ring == 'O'*/) {
            emptyEdist->Fill(mult);
-           
          }
-         else {
-           AliWarning("Something is wrong - wrong trigger");
-           continue;
+         //if(emptyCside && det == 3 && ring == 'I') {
+         //  emptyEdist->Fill(mult);
+         //}
+         
+         
+
+         if((emptyAside || emptyCside) && !(det == 3 /*&& ring == 'O'*/)/* && !(det == 2 && ring == 'O')*/) {
+           emptyEdist->Fill(mult);
          }
+         //if(emptyCside && det == 2 && ring == 'O') {
+         //  emptyEdist->Fill(mult);
+         // }
+        
+                 
+         //else {
+         //  AliWarning("Something is wrong - wrong trigger");
+         //  continue;
+         //}
         
          
        }
@@ -234,7 +256,7 @@ void AliFMDAnalysisTaskCollector::UserExec(Option_t */*option*/)
 //____________________________________________________________________
   
 void AliFMDAnalysisTaskCollector::Terminate(Option_t */*option*/) {
-  std::cout<<"Analysed "<<fEvents<<" events and "<<fEmptyEvents<<" empty"<<std::endl;
+  std::cout<<"Analysed "<<fEvents<<" events and "<<fEmptyEventsAside<<" empty from A side and "<<fEmptyEventsCside<<"  on the C side"<<std::endl;
   AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
  
     
@@ -245,8 +267,16 @@ void AliFMDAnalysisTaskCollector::Terminate(Option_t */*option*/) {
        Char_t ringChar = (ring == 0 ? 'I' : 'O');
        TH1F* hRingEdist = (TH1F*)fOutputList->FindObject(Form("ringFMD%d%c",det,ringChar));
        TH1F* hEmptyEdist = (TH1F*)fOutputList->FindObject(Form("emptyFMD%d%c",det,ringChar));
-       if(fEmptyEvents)
-         hEmptyEdist->Scale(1./(Float_t)fEmptyEvents);
+       if(fEmptyEventsAside && det == 3 /*&& ringChar == 'O'*/)
+         hEmptyEdist->Scale(1./(Float_t)fEmptyEventsAside);
+       
+       //if(fEmptyEventsAside && det == 3 && ringChar == 'I')
+       //  hEmptyEdist->Scale(1./(Float_t)fEmptyEventsCside);
+       //if(fEmptyEventsAside && det == 2 && ringChar == 'O') //A side
+       //  hEmptyEdist->Scale(1./(Float_t)fEmptyEventsAside);
+       
+       if((fEmptyEventsAside != 0 || fEmptyEventsCside !=0) && !(det == 3/* && ringChar == 'O'*/)/* && !(det == 2 && ringChar == 'O')*/)
+         hEmptyEdist->Scale(1./((Float_t)fEmptyEventsAside + ((Float_t)fEmptyEventsCside)));
        
        if(fEvents)
          hRingEdist->Scale(1./(Float_t)fEvents);
index d5461e0..15ce92d 100644 (file)
@@ -46,7 +46,8 @@ private:
   TObjArray*    fArray;
   TH1F*         fZvtxDist;
   Int_t         fEvents;
-  Int_t         fEmptyEvents;
+  Int_t         fEmptyEventsAside;
+  Int_t         fEmptyEventsCside;
   ClassDef(AliFMDAnalysisTaskCollector, 0); // Analysis task for FMD analysis
 };
  
index 4e5fe80..4534550 100644 (file)
@@ -17,6 +17,7 @@
 #include "AliAODHandler.h"
 #include "AliMCEventHandler.h"
 #include "AliStack.h"
+#include "AliFMDAnaCalibEnergyDistribution.h"
 #include "AliESDVertex.h"
 #include "TMath.h"
 #include "AliFMDAnaParameters.h"
@@ -65,6 +66,8 @@ void AliFMDAnalysisTaskDensity::CreateOutputObjects()
 {
   AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
   
+  
+  
   if(!fOutputList)
     fOutputList = new TList();
   fOutputList->SetName("density_list");
@@ -112,10 +115,12 @@ void AliFMDAnalysisTaskDensity::ConnectInputData(Option_t */*option*/)
     fVertex = (AliESDVertex*)GetInputData(1);
   }
 }
+//_____________________________________________________________________
+void AliFMDAnalysisTaskDensity::Init() {
+  //TFile f("/home/canute/ALICE/FMDanalysis/FirstAnalysis/energydistributions_0_0_1_0_0_0.root");
+  //fEnergyDistribution = dynamic_cast<AliFMDAnaCalibEnergyDistribution*>(f.Get("energydistributions"));
 
-
-
-
+}
 //_____________________________________________________________________
 void AliFMDAnalysisTaskDensity::Exec(Option_t */*option*/)
 {
@@ -124,9 +129,11 @@ void AliFMDAnalysisTaskDensity::Exec(Option_t */*option*/)
   
   //AliESDFMD*   fmd = fESD->GetFMDData();
   
+  
   Double_t vertex[3];
   fVertex->GetXYZ(vertex);
   // Z Vtx cut
+  
   if( TMath::Abs(vertex[2]) > pars->GetVtxCutZ()) {
     fStatus = kFALSE;
     return;
@@ -161,7 +168,12 @@ void AliFMDAnalysisTaskDensity::Exec(Option_t */*option*/)
      
       UShort_t nsec = (ir == 0 ? 20  : 40);
       UShort_t nstr = (ir == 0 ? 512 : 256);
-      
+      /*
+     TH1F* hEnergyDist       = pars->GetEmptyEnergyDistribution(det,ring);
+      TF1*  fitFunc           = hEnergyDist->GetFunction("FMDfitFunc");
+      TH1F* hSignalDist    = pars->GetRingEnergyDistribution(det, ring);
+      TF1*  fitFuncSignal  = hSignalDist->GetFunction("FMDfitFunc");
+      */
       for(UShort_t sec =0; sec < nsec;  sec++)  {
        for(UShort_t strip = 0; strip < nstr; strip++) {
          Float_t mult = fESD->Multiplicity(det,ring,sec,strip);
@@ -171,7 +183,7 @@ void AliFMDAnalysisTaskDensity::Exec(Option_t */*option*/)
          Float_t phi = pars->GetPhiFromSector(det,ring,sec);
          Float_t eta = pars->GetEtaFromStrip(det,ring,sec,strip,vertex[2]);
          
-         Float_t mult_cut = 0.25;//0.15;//m-2*s;//0.15;//0.2;//m-3*s;// 0.2;//0.01;//m-2*s;//0.2;
+         Float_t mult_cut = 0.3;//pars->GetMPV(det,ring,eta);// - pars->GetSigma(det,ring,eta);//0.25;//0.15;//m-2*s;//0.15;//0.2;//m-3*s;// 0.2;//0.01;//m-2*s;//0.2;
          // if(ring == 'I')
          //  mult_cut = 0.10;
          //Float_t mult_cut = pars->GetMPV(det,ring,eta) - 5*pars->GetSigma(det,ring,eta);
@@ -180,8 +192,7 @@ void AliFMDAnalysisTaskDensity::Exec(Option_t */*option*/)
            //proton + proton
            
            if(mult > mult_cut) {
-             nParticles = 1; 
-           
+             nParticles = 1.;
            }
          }
          else {
@@ -224,6 +235,27 @@ void AliFMDAnalysisTaskDensity::Exec(Option_t */*option*/)
            
          }
          
+         //std::cout<<det<<"   "<<ring<<"    "<<sec<<"    "<<strip<<"    "<<vertex[2]<<"   "<<eta<<std::endl;
+         
+         //Float_t signal = hSignalDist->Integral(hSignalDist->FindBin(0.5),hSignalDist->FindBin(2)) ;//pars->GetConstant(det,ring,eta);
+         //Float_t bkg =  hEnergyDist->GetBinContent(hEnergyDist->FindBin(mult));
+         //Float_t signal =  hSignalDist->GetBinContent(hSignalDist->FindBin(mult));
+         /*
+         if(fitFunc && fitFuncSignal && pars->IsRealData()) {
+           Float_t bkg = fitFunc->Eval(mult);
+           Float_t signal = fitFuncSignal->Eval(mult);
+           
+           if(signal > 0) {
+             correction = correction/(1-(bkg/signal));
+             //test
+             // if(TMath::Abs(eta)<3) correction = 1.15*correction;
+             // if(det == 2 && ring == 'I') correction = correction/(1-bkg/signal);
+             //if(det == 2 && ring == 'O') correction = correction/(1-bkg/signal);
+             //if(det == 3 && ring == 'I') correction = correction/(1-bkg/signal);
+             //if(det == 3 && ring == 'O') correction = correction/(1-bkg/signal);
+             
+           }
+           }*/
          if(correction) nParticles = nParticles / correction;
          if(nParticles > 0)
            hMult->Fill(eta,phi,nParticles);
@@ -238,6 +270,7 @@ void AliFMDAnalysisTaskDensity::Exec(Option_t */*option*/)
   
   }
   
+  
 
   if(fStandalone) {
     PostData(0, fOutputList); 
index 7b1f4f4..a8e97f6 100644 (file)
@@ -16,7 +16,7 @@ class AliESDEvent;
 class TChain;
 class AliAODEvent;
 class TF1;
-
+class AliFMDAnaCalibEnergyDistribution;
 
 /**
  * @ingroup FMD_ana
@@ -39,7 +39,7 @@ class AliFMDAnalysisTaskDensity : public AliAnalysisTask
     // Implementation of interface methods
     virtual void ConnectInputData(Option_t *option);
     virtual void CreateOutputObjects();
-    virtual void Init() {}
+  virtual void Init();
     virtual void LocalInit() {Init();}
     virtual void Exec(Option_t *option);
     virtual void Terminate(Option_t */*option*/) {}
@@ -53,17 +53,18 @@ class AliFMDAnalysisTaskDensity : public AliAnalysisTask
     Float_t GetPhiFromSector(UShort_t det, Char_t ring, UShort_t sec);
  private:
     
-    Int_t         fDebug;        //  Debug flag
-    TList*        fOutputList;
-    AliESDFMD*    fESD;
-    TObjString    fVertexString;
-    AliESDVertex* fVertex;
-    Bool_t        fStandalone;
-    Bool_t        fStatus;
+  Int_t         fDebug;        //  Debug flag
+  TList*        fOutputList;
+  AliESDFMD*    fESD;
+  TObjString    fVertexString;
+  AliESDVertex* fVertex;
+  Bool_t        fStandalone;
+  Bool_t        fStatus;
+  //AliFMDAnaCalibEnergyDistribution* fEnergyDistribution;
     // TF1*          fFuncPos;
     // TF1*          fFuncNeg;
     
-    ClassDef(AliFMDAnalysisTaskDensity, 0); // Analysis task for FMD analysis
+  ClassDef(AliFMDAnalysisTaskDensity, 0); // Analysis task for FMD analysis
 };
  
 #endif
index a554201..5bd7813 100644 (file)
@@ -51,7 +51,8 @@ void AliFMDAnalysisTaskSE::UserCreateOutputObjects()
   fSharing.SetFMDData(fmd);
   fSharing.SetVertex(vertex);
   fSharing.SetOutputList(fListOfHistos);
-
+  
+  fDensity.Init();
   fDensity.SetOutputList(densitylist);
   fDensity.SetInputESDFMD(fmd) ;
   fDensity.SetInputVertex(vertex);
@@ -94,9 +95,9 @@ void AliFMDAnalysisTaskSE::UserExec(Option_t */*option*/)
       fDndeta.Exec("");
       
     }
+    else return;
   }
-  else
-    return;
+  else return;
   
   //fListOfHistos = fBackground.GetOutputList();
   
@@ -105,9 +106,21 @@ void AliFMDAnalysisTaskSE::UserExec(Option_t */*option*/)
 //_____________________________________________________________________
 void AliFMDAnalysisTaskSE::Terminate(Option_t */*option*/)
 {
+  
+  TList* outputList = (TList*)GetOutputData(1);
+  
+  
+  fSharing.SetOutputList(outputList);
+  fBackground.SetHitList(outputList);
+  fDndeta.SetOutputList(outputList); 
+  fSharing.Terminate("");
   fBackground.Terminate("");
   fDndeta.Terminate("");
-
+  
+  
+  // TFile file("fmd_ana_histos_tmp.root","RECREATE");
+  //  fListOfHistos->Write();
+  // file.Close();
 }
 
 //_____________________________________________________________________
index b3b1c66..a2d5030 100644 (file)
@@ -51,12 +51,12 @@ public:                             //
   void         Print(Option_t* option="") const;
 private:
   
-  TList* fListOfHistos;
-  AliFMDAnalysisTaskSharing              fSharing;     // Sharing task
-  AliFMDAnalysisTaskDensity              fDensity;     // Density task
-  AliFMDAnalysisTaskBackgroundCorrection fBackground;  // Background task
-  AliFMDAnalysisTaskDndeta               fDndeta;      // dN/deta task
-  AliFMDAnaParameters*                   fParams;      // Analysis parameters
+  TList*                                 fListOfHistos; // Output list
+  AliFMDAnalysisTaskSharing              fSharing;      // Sharing task
+  AliFMDAnalysisTaskDensity              fDensity;      // Density task
+  AliFMDAnalysisTaskBackgroundCorrection fBackground;   // Background task
+  AliFMDAnalysisTaskDndeta               fDndeta;       // dN/deta task
+  AliFMDAnaParameters*                   fParams;       // Analysis parameters
 
   
   ClassDef(AliFMDAnalysisTaskSE, 1);
index 3a899bc..e8862ac 100644 (file)
 //#include "/home/canute/ALICE/AliRoot/PWG0/AliPWG0Helper.h"
 //#include "AliFMDParameters.h"
 #include "AliGenEventHeader.h"
+#include "AliGenPythiaEventHeader.h"
 #include "AliHeader.h"
 #include "AliStack.h"
 #include "AliMCParticle.h"
 #include "AliFMDStripIndex.h"
+#include "AliESDVZERO.h"
+#include "AliESDtrack.h"
 
 // This is the task to do the FMD sharing or hit merging.
 // It reads the input ESDFMD data and posts an ESDFMD object to
@@ -69,7 +72,8 @@ AliFMDAnalysisTaskSharing::AliFMDAnalysisTaskSharing(const char* name, Bool_t SE
     fStandalone(kTRUE),
     fEsdVertex(0),
     fStatus(kTRUE),
-    fLastTrackByStrip(0)
+    fLastTrackByStrip(0),
+    fEtotal(0)
 {
   // named constructor
   fStandalone = SE;
@@ -104,8 +108,12 @@ void AliFMDAnalysisTaskSharing::CreateOutputObjects()
                            hBg->GetXaxis()->GetXmax());
   hPrimary->Sumw2();
   fDiagList->Add(hPrimary);
+  TH1F* hXvtx = new TH1F("hXvtx","x vertex distribution",100,-2,2);
+  TH1F* hYvtx = new TH1F("hYvtx","y vertex distribution",100,-2,2);
   TH1F* hZvtx = new TH1F("hZvtx","z vertex distribution",pars->GetNvtxBins(),-1*pars->GetVtxCutZ(),pars->GetVtxCutZ());
   
+  fDiagList->Add(hXvtx);
+  fDiagList->Add(hYvtx);
   fDiagList->Add(hZvtx);
   
   TH1F* hPrimVertexBin = 0;
@@ -154,10 +162,87 @@ void AliFMDAnalysisTaskSharing::CreateOutputObjects()
       
     }
   }
+  TH2F*  hCorrelationFMDSPDhits = new TH2F("hCorrelationFMDSPDhits","hCorrelationFMDSPDhits;SPD;FMD ",100,0,200,100,0,500);
+  TH2F*  hCorrelationFMDSPD = new TH2F("hCorrelationFMDSPD","hCorrelationFMDSPD;SPD;FMD ",100,0,200,100,0,500);
+  TH2F*  hCorrelationFMD1SPD = new TH2F("hCorrelationFMD1SPD","hCorrelationFMD1SPD;SPD;FMD1 ",100,0,200,100,0,200);
+  TH2F*  hCorrelationFMD2ISPD = new TH2F("hCorrelationFMD2ISPD","hCorrelationFMD2ISPD;SPD;FMD2I ",100,0,200,100,0,200);
+  TH2F*  hCorrelationFMD2OSPD = new TH2F("hCorrelationFMD2OSPD","hCorrelationFMD2OSPD;SPD;FMD2O ",100,0,200,100,0,200);
+  TH2F*  hCorrelationFMD3ISPD = new TH2F("hCorrelationFMD3ISPD","hCorrelationFMD3ISPD;SPD;FMD3I ",100,0,200,100,0,200);
+  TH2F*  hCorrelationFMD3OSPD = new TH2F("hCorrelationFMD3OSPD","hCorrelationFMD3OSPD;SPD;FMD3O ",100,0,200,100,0,200);
+  TH2F*  hCorrelationFMDGoodtracks = new TH2F("hCorrelationFMDGoodtracks","hCorrelationGoodtracks;good tracks;FMD ",100,0,200,100,0,200);
+  TH2F*  hCorrelationFMDBadtracks = new TH2F("hCorrelationFMDBadtracks","hCorrelationBadtracks;bad tracks;FMD ",100,0,200,100,0,200);
+  TH2F*  hCorrelationGoodbadtracks = new TH2F("hCorrelationGoodbadtracks","hCorrelationGoodbadtracks;good tracks;bad tracks ",100,0,200,100,0,200);
+  TH2F*  hCorrelationSPDTracklets = new TH2F("hCorrelationSPDTracklets","hCorrelationSPDTracklets;hits ; tracklets ",100,0,500,100,0,200);
+  TH2F*  hCorrelationClustersTracklets = new TH2F("hCorrelationClustersTracklets","hCorrelationClustersTracklets;clusters ; tracklets ",500,0,500,100,0,100);
+  TH2F*  hCorrelationHitsRadius = new TH2F("hCorrelationHitsRadius","hCorrelationHitsRadius;hits ; radius ",100,0,500,100,0,10);
+  TH2F*  hCorrelationHitsX = new TH2F("hCorrelationHitsX","hCorrelationHitsX;hits ; X ",100,0,500,100,-5,5);
+  TH2F*  hCorrelationHitsY = new TH2F("hCorrelationHitsY","hCorrelationHitsY;hits ; Y ",100,0,500,100,-5,5);
+  fDiagList->Add(hCorrelationHitsRadius);
+  fDiagList->Add(hCorrelationHitsX);
+  fDiagList->Add(hCorrelationHitsY);
+  fDiagList->Add(hCorrelationFMDSPD);
+  fDiagList->Add(hCorrelationFMD1SPD);
+  fDiagList->Add(hCorrelationFMD2ISPD);
+  fDiagList->Add(hCorrelationFMD2OSPD);
+  fDiagList->Add(hCorrelationFMD3ISPD);
+  fDiagList->Add(hCorrelationFMD3OSPD);
+  fDiagList->Add(hCorrelationFMDGoodtracks);
+  fDiagList->Add(hCorrelationFMDBadtracks);
+  fDiagList->Add(hCorrelationGoodbadtracks);
+fDiagList->Add(hCorrelationFMDSPDhits);
+  fDiagList->Add(hCorrelationClustersTracklets);
+  fDiagList->Add(hCorrelationSPDTracklets);
+  TH2F*  hCorrelationFMD23 = new TH2F("hCorrelationFMD23","hCorrelationFMD23;FMD2 ;FMD3 ",100,0,500,100,0,500);
+  TH2F*  hCorrelationFMD2diff23 = new TH2F("hCorrelationFMD2diff23","hCorrelationFMD2diff23;FMD2 ;diff FMD23 ",100,0,100,100,0,100);
+  TH2F*  hCorrelationFMD3diff23 = new TH2F("hCorrelationFMD3diff23","hCorrelationFMD3diff23;FMD3 ;diff FMD23 ",100,0,100,100,0,100);
+  TH2F*  hCorrelationFMD1diff13 = new TH2F("hCorrelationFMD1diff13","hCorrelationFMD1diff13;FMD1 ;diff FMD13 ",100,0,100,100,0,100);
+  TH2F*  hCorrelationFMD1diff12 = new TH2F("hCorrelationFMD1diff12","hCorrelationFMD1diff12;FMD1 ;diff FMD12 ",100,0,100,100,0,100);
+  TH2F*  hCorrelationFMD12 = new TH2F("hCorrelationFMD12","hCorrelationFMD12;FMD1 ;FMD2 ",100,0,500,100,0,500);
+  TH2F* hCorrelationFMDBgCand = new TH2F("hCorrelationFMDBgCand","hCorrelationFMDBgCand;Bg Tr ;FMD ",100,0,100,500,0,500);
+  
+  TH2F* hCorrelationFMDFlatTr = new TH2F("hCorrelationFMDFlatTr","hCorrelationFMDFlatTr;Bg Tr ;FMD ",100,0,100,500,0,500);
+  TH2F* hCorrelationFMDRatioFlatTr = new TH2F("hCorrelationFMDRatioFlatTr","hCorrelationFMDRatioFlatTr;Bg Tr ;FMD ",100,0,1,500,0,500);
+  fDiagList->Add(hCorrelationFMDBgCand);
+  fDiagList->Add(hCorrelationFMDFlatTr);
+  fDiagList->Add(hCorrelationFMDRatioFlatTr);
+  TH2F* hCorrelationFMDBgCandRelative = new TH2F("hCorrelationFMDBgCandRelative","hCorrelationFMDBgCandRelative;Bg Tr ;FMD ",100,0,2,500,0,500);
+  fDiagList->Add(hCorrelationFMDBgCandRelative);
+  fDiagList->Add(hCorrelationFMD2diff23);
+  fDiagList->Add(hCorrelationFMD3diff23);
+  fDiagList->Add(hCorrelationFMD1diff13);
+  fDiagList->Add(hCorrelationFMD1diff12);
+  fDiagList->Add(hCorrelationFMD23);
+  fDiagList->Add(hCorrelationFMD12);
+  TH2F*  hTimeCorrelation = new TH2F("hCorrelationTime","hCorrelationTime ; time ; FMD hits",500,0,500,100,0,200);
+  fDiagList->Add(hTimeCorrelation);
+  TH1F*  hHitDistribution = new TH1F("hFMDHitDistribution","hFMDHitDistribution ; FMD hits",500,0,500);
+
+  TH1F*  hHitDistributionFMD1 = new TH1F("hFMDHitDistributionFMD1","hFMDHitDistributionFMD1 ; FMD1 hits",500,0,500);
+  TH1F*  hHitDistributionFMD2I = new TH1F("hFMDHitDistributionFMD2I","hFMDHitDistributionFMD2I ; FMD2I hits",500,0,500);
+  TH1F*  hHitDistributionFMD2O = new TH1F("hFMDHitDistributionFMD2O","hFMDHitDistributionFMD2O ; FMD2O hits",500,0,500);
+  TH1F*  hHitDistributionFMD3I = new TH1F("hFMDHitDistributionFMD3I","hFMDHitDistributionFMD3I ; FMD3I hits",500,0,500);
+  TH1F*  hHitDistributionFMD3O = new TH1F("hFMDHitDistributionFMD3O","hFMDHitDistributionFMD3O ; FMD3O hits",500,0,500);
+  TH1F*  hTrVtxDistribution = new TH1F("hTrVtxDistribution","hTrVtxDistribution ; TrVtx",200,-500,500);
+  TH1F*  hTrEtaDistribution = new TH1F("hTrEtaDistribution","hTrEtaDistribution ; TrEta",200,-9,9);
+  TH1F*  hTrEtaDistribution2 = new TH1F("hTrEtaDistribution2","hTrEtaDistribution2 ; TrEta",200,-9,9);
+  
+ TH1F*  hFlatTracks = new TH1F("hFlatTracks","hFlatTracks ; Horizontal tracks",100,0,100);
+ fDiagList->Add(hTrVtxDistribution);
+ fDiagList->Add(hTrEtaDistribution);
+ fDiagList->Add(hTrEtaDistribution2);
+  fDiagList->Add(hFlatTracks);
+  fDiagList->Add(hHitDistribution);
+  fDiagList->Add(hHitDistributionFMD1);
+  fDiagList->Add(hHitDistributionFMD2I);
+  fDiagList->Add(hHitDistributionFMD2O);
+  fDiagList->Add(hHitDistributionFMD3I);
+  fDiagList->Add(hHitDistributionFMD3O);
   TH1F*  nMCevents = new TH1F("nMCEventsNoCuts","nMCEventsNoCuts",pars->GetNvtxBins(),0,pars->GetNvtxBins());
   
   fDiagList->Add(nMCevents);
-  
+    
+    
 }
 //_____________________________________________________________________
 void AliFMDAnalysisTaskSharing::ConnectInputData(Option_t */*option*/)
@@ -177,6 +262,10 @@ void AliFMDAnalysisTaskSharing::Exec(Option_t */*option*/)
   
   foutputESDFMD->Clear();
   
+  Int_t delta = fESD->GetOrbitNumber() - fLastOrbit;
+  fLastOrbit = fESD->GetOrbitNumber();
+  
+  
   AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
   Double_t vertex[3];
   Bool_t vtxStatus = pars->GetVertex(fESD,vertex);
@@ -185,37 +274,45 @@ void AliFMDAnalysisTaskSharing::Exec(Option_t */*option*/)
   // Process primaries here to get true MC distribution
   if(pars->GetProcessPrimary())
     ProcessPrimary();
+  const AliMultiplicity* testmult = fESD->GetMultiplicity();
+  Int_t nTrackLets = testmult->GetNumberOfTracklets();
+  TH2F*  hCorrelationClustersTracklets = (TH2F*)fDiagList->FindObject("hCorrelationClustersTracklets");
+  hCorrelationClustersTracklets->Fill(testmult->GetNumberOfSingleClusters(),nTrackLets);  
   
-  Bool_t isTriggered = pars->IsEventTriggered(fESD);
   
-  if(!isTriggered) {
-    fStatus = kFALSE;
-    return;
-   }
-   else
-     fStatus = kTRUE;
   
-  if(!vtxStatus) {
+  Bool_t isTriggered = pars->IsEventTriggered(fESD);
+  
+  if(!isTriggered || !vtxStatus ) {
     fStatus = kFALSE;
     return;
   }
   else
     fStatus = kTRUE;
   
+  TH1F* hXvtx = (TH1F*)fDiagList->FindObject("hXvtx");
+  if(vertex[0] != 0) hXvtx->Fill(vertex[0]);
+  TH1F* hYvtx = (TH1F*)fDiagList->FindObject("hYvtx");
+  if(vertex[1] != 0) hYvtx->Fill(vertex[1]);
   TH1F* hZvtx = (TH1F*)fDiagList->FindObject("hZvtx");
   hZvtx->Fill(vertex[2]);
-  
-  const AliMultiplicity* testmult = fESD->GetMultiplicity();
+  //const AliMultiplicity* testmult = fESD->GetMultiplicity();
+  //std::cout<<vertex[2]<<std::endl;
+  //Int_t nTrackLets = testmult->GetNumberOfTracklets();
   
-  Int_t nTrackLets = testmult->GetNumberOfTracklets();
+  if( TMath::Abs(vertex[2]) > pars->GetVtxCutZ()) {
+    fStatus = kFALSE;
+    return;
+  }
+    
   if(nTrackLets < 1000) foutputESDFMD->SetUniqueID(kTRUE);
   else foutputESDFMD->SetUniqueID(kFALSE);
   
   AliESDFMD* fmd = fESD->GetFMDData();
   
   if (!fmd) return;
-  Int_t nHits = 0;
+  Int_t nHits[3][2] = {{0,0},{0,0},{0,0}};
+
   for(UShort_t det=1;det<=3;det++) {
     Int_t nRings = (det==1 ? 1 : 2);
     for (UShort_t ir = 0; ir < nRings; ir++) {
@@ -235,11 +332,10 @@ void AliFMDAnalysisTaskSharing::Exec(Option_t */*option*/)
          
          if(mult == AliESDFMD::kInvalidMult || mult == 0) continue;
          
-         //Double_t eta  = EtaFromStrip(det,ring,sec,strip,vertex[2]);//fmd->Eta(det,ring,sec,strip);
+        
          //Double_t eta = fmd->Eta(det,ring,sec,strip);
          Float_t eta = pars->GetEtaFromStrip(det,ring,sec,strip,vertex[2]);
-         //std::cout<<EtaFromStrip(det,ring,sec,strip,vertex[2]) <<"    "<<fmd->Eta(det,ring,sec,strip)<<std::endl;
-         
+                 
          hEdist->Fill(mult);
          if(fmd->IsAngleCorrected())
            mult = mult/TMath::Cos(Eta2Theta(eta));
@@ -259,9 +355,11 @@ void AliFMDAnalysisTaskSharing::Exec(Option_t */*option*/)
            }
          
          Float_t mergedEnergy = GetMultiplicityOfStrip(mult,eta,prevE,nextE,det,ring,sec,strip);
-
-         if(mergedEnergy > 0 )
-           nHits++;
+         //if(mult> (pars->GetMPV(det,ring,eta) - pars->GetSigma(det,ring,eta))) 
+         //  mergedEnergy = mult;
+         //else mergedEnergy = 0;
+         if(mergedEnergy > 0.3 )
+           nHits[det-1][ir]++;
          foutputESDFMD->SetMultiplicity(det,ring,sec,strip,mergedEnergy);
          foutputESDFMD->SetEta(det,ring,sec,strip,eta);
          
@@ -269,8 +367,153 @@ void AliFMDAnalysisTaskSharing::Exec(Option_t */*option*/)
       }
     }
   }
+  //cluster cut
+  //if(testmult->GetNumberOfSingleClusters() > 15 + nTrackLets)
+  //  {fStatus = kFALSE; std::cout<<"FMD : "<<nHits[0][0]<<"  "<<nHits[1][0]<<"   "<<nHits[1][1]<<"   "<<nHits[2][0]<<"   "<<nHits[2][1]<<" tracks  "<<testmult->GetNumberOfSingleClusters()<<"   "<<nTrackLets<<std::endl; return;}
+
+  TH2F*  hCorrelationFMD23  = (TH2F*)fDiagList->FindObject("hCorrelationFMD23");
+  TH2F*  hCorrelationFMD12  = (TH2F*)fDiagList->FindObject("hCorrelationFMD12");
+  TH2F*  hCorrelationFMD2diff23  = (TH2F*)fDiagList->FindObject("hCorrelationFMD2diff23");
+  TH2F*  hCorrelationFMD3diff23  = (TH2F*)fDiagList->FindObject("hCorrelationFMD3diff23");
+  TH2F*  hCorrelationFMD1diff13  = (TH2F*)fDiagList->FindObject("hCorrelationFMD1diff13");  
+  TH2F*  hCorrelationFMD1diff12  = (TH2F*)fDiagList->FindObject("hCorrelationFMD1diff12");
+
+  TH2F*  hCorrelationFMDSPD = (TH2F*)fDiagList->FindObject("hCorrelationFMDSPD");
+  TH2F*  hCorrelationFMD1SPD = (TH2F*)fDiagList->FindObject("hCorrelationFMD1SPD");
+  TH2F*  hCorrelationFMD2ISPD = (TH2F*)fDiagList->FindObject("hCorrelationFMD2ISPD");
+  TH2F*  hCorrelationFMD2OSPD = (TH2F*)fDiagList->FindObject("hCorrelationFMD2OSPD");
+  TH2F*  hCorrelationFMD3ISPD = (TH2F*)fDiagList->FindObject("hCorrelationFMD3ISPD");
+  TH2F*  hCorrelationFMD3OSPD = (TH2F*)fDiagList->FindObject("hCorrelationFMD3OSPD");
+  
+  TH2F*  hCorrelationFMDSPDhits = (TH2F*)fDiagList->FindObject("hCorrelationFMDSPDhits");
+
+  TH2F*  hCorrelationSPDTracklets = (TH2F*)fDiagList->FindObject("hCorrelationSPDTracklets");
+  TH2F*  hTimeCorrelation   = (TH2F*)fDiagList->FindObject("hCorrelationTime");
+  TH2F*  hHitsRadius   = (TH2F*)fDiagList->FindObject("hCorrelationHitsRadius");
+  TH2F*  hHitsX   = (TH2F*)fDiagList->FindObject("hCorrelationHitsX");
+  TH2F*  hHitsY   = (TH2F*)fDiagList->FindObject("hCorrelationHitsY");
+  TH1F*  hFMDHitDistribution   = (TH1F*)fDiagList->FindObject("hFMDHitDistribution");
+  TH1F*  hFMD1HitDistribution   = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD1");  
+  TH1F*  hFMD2IHitDistribution   = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD2I");
+  TH1F*  hFMD2OHitDistribution   = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD2O");
+  TH1F*  hFMD3IHitDistribution   = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD3I");
+  TH1F*  hFMD3OHitDistribution   = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD3O"); 
+  TH1F*  hCorrelationFMDGoodtracks = (TH1F*)fDiagList->FindObject("hCorrelationFMDGoodtracks");
+  TH1F*  hCorrelationFMDBadtracks = (TH1F*)fDiagList->FindObject("hCorrelationFMDBadtracks");
+  TH1F*  hCorrelationGoodbadtracks = (TH1F*)fDiagList->FindObject("hCorrelationGoodbadtracks");
+   TH2F*  hCorrelationFMDBgCand = (TH2F*)fDiagList->FindObject("hCorrelationFMDBgCand");
+   TH2F*  hCorrelationFMDBgCandRelative = (TH2F*)fDiagList->FindObject("hCorrelationFMDBgCandRelative");
+   
+   TH2F* hCorrelationFMDFlatTr = (TH2F*)fDiagList->FindObject("hCorrelationFMDFlatTr");
+   TH2F* hCorrelationFMDRatioFlatTr = (TH2F*)fDiagList->FindObject("hCorrelationFMDRatioFlatTr");
+   TH1F* hTrVtxDistribution = (TH1F*)fDiagList->FindObject("hTrVtxDistribution");
+   TH1F* hTrEtaDistribution = (TH1F*)fDiagList->FindObject("hTrEtaDistribution");
+   TH1F* hTrEtaDistribution2 = (TH1F*)fDiagList->FindObject("hTrEtaDistribution2");
+  hCorrelationFMDSPD->Fill(nTrackLets,nHits[0][0]+nHits[1][0]+nHits[1][1]+nHits[2][0]+nHits[2][1]);
+  TH1F* hFlatTracks = (TH1F*)fDiagList->FindObject("hFlatTracks");
+  hCorrelationFMD1SPD->Fill(nTrackLets,nHits[0][0]);
+  hCorrelationFMD2ISPD->Fill(nTrackLets,nHits[1][0]);
+  hCorrelationFMD2OSPD->Fill(nTrackLets,nHits[1][1]);
+  hCorrelationFMD3ISPD->Fill(nTrackLets,nHits[2][0]);
+  hCorrelationFMD3OSPD->Fill(nTrackLets,nHits[2][1]);
+  hCorrelationFMDSPDhits->Fill(testmult->GetNumberOfFiredChips(0),nHits[0][0]+nHits[1][0]+nHits[1][1]+nHits[2][0]+nHits[2][1]);
+  hCorrelationSPDTracklets->Fill(testmult->GetNumberOfFiredChips(0),nTrackLets);
+  
+  hTimeCorrelation->Fill(delta,nHits[0][0]+nHits[1][0]+nHits[1][1]+nHits[2][0]+nHits[2][1]);
+  hCorrelationFMD23->Fill(nHits[1][0]+nHits[1][1],nHits[2][0]+nHits[2][1]);
+  hCorrelationFMD12->Fill(nHits[0][0],nHits[1][0]+nHits[1][1]);
+  
+  //  if(TMath::Abs(nHits[1]-nHits[2]) > 15 && (nHits[1]+nHits[2]) > 35) {fStatus = kFALSE; std::cout<<"difference : "<<TMath::Abs(nHits[1]-nHits[2])<<std::endl; return;}
+  
+  // if(testmult->GetNumberOfFiredChips(0))
+  //  if(testmult->GetNumberOfFiredChips(0) > 15 && ((Float_t)nTrackLets / (Float_t)testmult->GetNumberOfFiredChips(0)) < 0.4)
+  //   {fStatus = kFALSE; std::cout<<nTrackLets<<"   "<<testmult->GetNumberOfFiredChips(0)<<"   "<<nHits[0]<<"   "<<nHits[1]<<"   "<<nHits[2]<<std::endl; return;}
+  
+  
+  Float_t diff23 = (Float_t)TMath::Abs(nHits[2][1] + nHits[2][0] - nHits[1][1] - nHits[1][0]);
+  Float_t diff13 = TMath::Abs(nHits[2][1] + nHits[2][0] - nHits[1][1] - nHits[1][0] - nHits[0][0]);
+  Float_t diff12 = TMath::Abs(nHits[1][0] - nHits[0][0]);
+  
+  hCorrelationFMD1diff12->Fill(nHits[0][0], diff12);
+  hCorrelationFMD1diff13->Fill(nHits[0][0], diff13);
+  hCorrelationFMD2diff23->Fill(nHits[1][1], diff23);
+  hCorrelationFMD3diff23->Fill(nHits[2][1], diff23);
+  
+  //
+  Float_t nTotalFMDhits = nHits[0][0]+nHits[1][0]+nHits[1][1]+nHits[2][0]+nHits[2][1] ;
+   Float_t radius = TMath::Sqrt(TMath::Power(vertex[0] + 0.03715,2) + TMath::Power(vertex[1] - 0.1659,2));
+  
+  if(vertex[1] !=0 || vertex[1] !=0) {
+    hHitsRadius->Fill(nTotalFMDhits,radius);
+    hHitsX->Fill(nTotalFMDhits,vertex[0]);
+    hHitsY->Fill(nTotalFMDhits,vertex[1]); }
   
+  hFMDHitDistribution->Fill(nTotalFMDhits);
+  hFMD1HitDistribution->Fill(nHits[0][0]);
+  hFMD2IHitDistribution->Fill(nHits[1][0]);
+  hFMD2OHitDistribution->Fill(nHits[1][1]);
+  hFMD3IHitDistribution->Fill(nHits[2][0]);
+  hFMD3OHitDistribution->Fill(nHits[2][0]);
+  
+  // if(radius > 0.5) {fStatus = kFALSE; std::cout<<"FMD : "<<nTotalFMDhits<<std::endl; foutputESDFMD->Clear(); return;}
+  //if(TMath::Abs(vertex[1] - 0.1659) > 0.1 ) {fStatus = kFALSE; std::cout<<"FMD : "<<nTotalFMDhits<<std::endl; foutputESDFMD->Clear(); return;}
+  
+   if(nTrackLets < pars->GetLowSPDLimit() || nTrackLets > pars->GetHighSPDLimit())
+    {fStatus = kFALSE; std::cout<<nTrackLets<<"   "<<"   "<<nHits[0][0]<<"  "<<nHits[1][0]<<"   "<<nHits[1][1]<<"   "<<nHits[2][0]<<"   "<<nHits[2][1]<<std::endl; return;}
+  
+    
+   AliESDtrack* track = 0;
+  Int_t ntracks = fESD->GetNumberOfTracks();
+  Float_t ngood =0, nbad = 0;
+  //std::cout<<" Primary vtx : "<<vertex[0]<<"   "<<vertex[1]<<"   "<<vertex[2]<<"  "<<nTotalFMDhits<<std::endl;
+  Int_t nBgCandidates = 0;
+  Float_t nFlat = 0;
+  for(Int_t i=0;i<ntracks;i++) {
+    track = fESD->GetTrack(i);
+    //std::cout<<track->GetX()-vertex[0]<<"   "<<track->GetY()-vertex[1]<<"   "<<track->GetZ()-vertex[2]<<std::endl;
+    //std::cout<<track->GetX()<<"   "<<track->GetY()<<"   "<<track->GetZ()<<std::endl;
+    hTrVtxDistribution->Fill(track->GetZ());
+    
+    if(TMath::Abs( track->GetZ()) > 50 &&  TMath::Abs(track->GetZ()) < 300) { // && TMath::Abs( track->GetY()) < 1)
+      nBgCandidates++;
+      hTrEtaDistribution->Fill(track->Eta()); 
+    }
+    
+    
+    
+    if(TMath::Abs(track->GetX()-vertex[0]) > 0.3 || TMath::Abs(track->GetY()-vertex[1]) > 0.3  || TMath::Abs(track->GetZ()-vertex[2]) > 0.3)  {
+      nbad++;
+      hTrEtaDistribution2->Fill(track->Eta()); }
+    else ngood++;
+    
+    if(TMath::Abs(track->Pt()) < 0.1)
+      nFlat++;
+  }
+  
+  Float_t ratioFlat = 0;
+  if(fESD->GetNumberOfTracks())
+    ratioFlat = nFlat/(Float_t)fESD->GetNumberOfTracks();
+  hCorrelationFMDFlatTr->Fill(nFlat,nTotalFMDhits);
+  hCorrelationFMDRatioFlatTr->Fill(ratioFlat,nTotalFMDhits);
+  hFlatTracks->Fill(nFlat);
    
+  // std::cout<<fESD->GetT0zVertex()<<"   "<<vertex[2]<<std::endl;
+  Float_t ratioBg = 0;
+  //if(fESD->GetNumberOfTracks() > 0)
+  
+  if(fESD->GetNumberOfTracks() > 0)
+    ratioBg = (Float_t)nBgCandidates/(Float_t)fESD->GetNumberOfTracks();
+  hCorrelationFMDBgCand->Fill(nBgCandidates,nTotalFMDhits);
+  hCorrelationFMDBgCandRelative->Fill(ratioBg,nTotalFMDhits);
+  
+  
+   Float_t ratio =  (nbad > 0 ? ngood / nbad  : 0);
+  
+  hCorrelationFMDGoodtracks->Fill(ngood,nTotalFMDhits);
+  hCorrelationFMDBadtracks->Fill(nbad,nTotalFMDhits);
+  hCorrelationGoodbadtracks->Fill(ngood,nbad);
+    
   if(fStandalone) {
     PostData(0, foutputESDFMD); 
     PostData(1, fEsdVertex); 
@@ -292,18 +535,28 @@ Float_t AliFMDAnalysisTaskSharing::GetMultiplicityOfStrip(Float_t mult,
  
   Float_t mergedEnergy = 0;
   //Float_t nParticles = 0;
-  Float_t cutLow  = 0.25;//0.15;
+  Float_t cutLow  = 0.3;//0.15;
+  
+  Float_t cutHigh = pars->GetMPV(det,ring,eta) - 2*pars->GetSigma(det,ring,eta);
   // if(ring == 'I')
   //  cutLow = 0.1;
   
   //cutLow = 0;
   //AliFMDParameters* recopars = AliFMDParameters::Instance();
   //cutLow = (5*recopars->GetPedestalWidth(det,ring,sec,strip))/(recopars->GetPulseGain(det,ring,sec,strip)*recopars->GetDACPerMIP());
+  //if(foutputESDFMD->GetUniqueID() == kFALSE ) {
+  
+  if(mult > 12 || mult < cutLow)
+    {
+      //   std::cout<<"rejecting hit in FMD "<<det<<" "<<ring<<std::endl;
+      fSharedThis      = kFALSE;
+      fSharedPrev      = kFALSE;
+      return 0;
+    }
+  
   
   
   
-  Float_t cutHigh = pars->GetMPV(det,ring,eta) - 3*pars->GetSigma(det,ring,eta);
-   
   // Float_t cutPart = pars->GetMPV(det,ring,eta) - 5*pars->GetSigma(det,ring,eta);
   Float_t totalE  = mult;
   
@@ -326,13 +579,7 @@ Float_t AliFMDAnalysisTaskSharing::GetMultiplicityOfStrip(Float_t mult,
       fSharedPrev      = kFALSE;
       return 0;
     }
-  if(mult > 15)
-    {
-      //   std::cout<<"rejecting hit in FMD "<<det<<" "<<ring<<std::endl;
-      fSharedThis      = kFALSE;
-      fSharedPrev      = kFALSE;
-      return 0;
-    }
+  
   
   if(prevE > cutLow && prevE < cutHigh && !fSharedPrev ) {
     totalE += prevE;
@@ -342,29 +589,90 @@ Float_t AliFMDAnalysisTaskSharing::GetMultiplicityOfStrip(Float_t mult,
     totalE += nextE;
     fSharedThis      = kTRUE;
   }
+  totalE = totalE*TMath::Cos(Eta2Theta(eta));
   TH1F* hEdist = (TH1F*)fDiagList->FindObject(Form("Edist_after_sharing_FMD%d%c",det,ring));
-  hEdist->Fill(totalE);
+  if(totalE > cutLow)
+    hEdist->Fill(totalE);
+  
   
-  totalE = totalE*TMath::Cos(Eta2Theta(eta));
   if(totalE > 0) {
     
     mergedEnergy = totalE;
     fSharedPrev      = kTRUE;
     // if(det == 1 && ring =='I')
-    // std::cout<<Form("Merged signals %f %f %f into %f , %f in strip %d, sec %d, ring %c, det %d",prevE, mult, nextE, totalE/TMath::Cos(Eta2Theta(eta)),totalE,strip,sec,ring,det )<<std::endl;
   }
-    else{// if(totalE > 0) {
-      //if(det == 3 && ring =='I')
-      //       std::cout<<Form("NO HIT  for  %f %f %f into %f , %f in strip %d, sec %d, ring %c, det %d, cuts %f , %f",prevE, mult, nextE, totalE/TMath::Cos(Eta2Theta(eta)),totalE,strip,sec,ring,det,cutPart,cutHigh )<<std::endl;
-    fSharedThis      = kFALSE;
-    fSharedPrev      = kFALSE;
+  else{
+      fSharedThis      = kFALSE;
+      fSharedPrev      = kFALSE;
   }
+
   // mergedEnergy = mult;
   
+
+  /*   }
+else {
+    TH1F* hEdist = (TH1F*)fDiagList->FindObject(Form("Edist_after_sharing_FMD%d%c",det,ring));
+    if(mult > cutLow)
+      fEtotal+=mult;
+    if(mult < cutLow) {
+      mergedEnergy = fEtotal;
+      fEtotal = 0;
+      hEdist->Fill(mergedEnergy);
+      
+    }
+    
+   }*/
+  
   return mergedEnergy; 
   //}  
 }
+//_____________________________________________________________________
+void AliFMDAnalysisTaskSharing::Terminate(Option_t* /* option*/) {
+  
+  TH1F*  hFMDHitDistribution   = (TH1F*)fDiagList->FindObject("hFMDHitDistribution");
+  TH1F*  hFMD1HitDistribution   = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD1");
+  TH1F*  hFMD2IHitDistribution   = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD2I");
+  TH1F*  hFMD2OHitDistribution   = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD2O");
+  TH1F*  hFMD3IHitDistribution   = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD3I");
+  TH1F*  hFMD3OHitDistribution   = (TH1F*)fDiagList->FindObject("hFMDHitDistributionFMD3O");
+  
+  TH1F* hZvtx = (TH1F*)fDiagList->FindObject("hZvtx");
+  for(UShort_t det=1;det<=3;det++) {
+    Int_t nRings = (det==1 ? 1 : 2);
+    for (UShort_t ir = 0; ir < nRings; ir++) {
+      Char_t   ring = (ir == 0 ? 'I' : 'O');
+      
+      TH1F* hEdist = (TH1F*)fDiagList->FindObject(Form("Edist_before_sharing_FMD%d%c",det,ring));
+      TH1F* hEdistAfter = (TH1F*)fDiagList->FindObject(Form("Edist_after_sharing_FMD%d%c",det,ring));
+      if(hZvtx->GetEntries()) {
+       hEdist->Scale(1./(Float_t)hZvtx->GetEntries());
+       hEdistAfter->Scale(1./(Float_t)hZvtx->GetEntries());
+      }
+      
+    }
+    
+  } 
+  TH1F* hFlatTracks = (TH1F*)fDiagList->FindObject("hFlatTracks");
+  TH1F* hTrVtxDistribution = (TH1F*)fDiagList->FindObject("hTrVtxDistribution");
+  TH1F* hTrEtaDistribution = (TH1F*)fDiagList->FindObject("hTrEtaDistribution");
+  TH1F* hTrEtaDistribution2 = (TH1F*)fDiagList->FindObject("hTrEtaDistribution2");
+  if(hZvtx->GetEntries()) {
+    hFMDHitDistribution->Scale(1./(Float_t)hZvtx->GetEntries());
+    hFMD1HitDistribution->Scale(1./(Float_t)hZvtx->GetEntries());
+    hFMD2IHitDistribution->Scale(1./(Float_t)hZvtx->GetEntries());
+    hFMD2OHitDistribution->Scale(1./(Float_t)hZvtx->GetEntries());
+    hFMD3IHitDistribution->Scale(1./(Float_t)hZvtx->GetEntries());
+    hFMD3OHitDistribution->Scale(1./(Float_t)hZvtx->GetEntries());
+    hFlatTracks->Scale(1./(Float_t)hZvtx->GetEntries());
+    hTrVtxDistribution->Scale(1./(Float_t)hZvtx->GetEntries());
+    hTrEtaDistribution->Scale(1./(Float_t)hZvtx->GetEntries());
+    hTrEtaDistribution2->Scale(1./(Float_t)hZvtx->GetEntries());
+  }
+  
 
+  
+  
+}
 //_____________________________________________________________________
 Float_t AliFMDAnalysisTaskSharing::Eta2Theta(Float_t eta) const{
   //convert the eta of a strip to a theta
@@ -400,7 +708,24 @@ void AliFMDAnalysisTaskSharing::ProcessPrimary() {
   AliHeader* header            = mcEvent->Header();
   AliGenEventHeader* genHeader = header->GenEventHeader();
   
+  AliGenPythiaEventHeader* pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
   
+  if (!pythiaGenHeader) {
+    std::cout<<" no pythia header!"<<std::endl;
+    return; 
+       }
+
+       
+  Int_t pythiaType = pythiaGenHeader->ProcessType();
+  
+  //if(pythiaType==92||pythiaType==93){
+    //  std::cout<<"single diffractive"<<std::endl;
+    //  return;
+    // }
+  // if(pythiaType==94){
+  //  std::cout<<"double diffractive"<<std::endl;
+  //  return;
+  // }
   
   TArrayF vertex;
   genHeader->PrimaryVertex(vertex);
index f5193db..49d8313 100644 (file)
@@ -47,7 +47,7 @@ class AliFMDAnalysisTaskSharing : public AliAnalysisTask
     virtual void Init() {}
     virtual void LocalInit() {Init();}
     virtual void Exec(Option_t */*option*/);
-    virtual void Terminate(Option_t* /* option*/) {}
+    virtual void Terminate(Option_t* /* option*/);
     virtual void SetDebugLevel(Int_t level) {fDebug = level;}
     Float_t GetMultiplicityOfStrip(Float_t mult, Float_t eta, Float_t Eprev, Float_t Enext, UShort_t   det, Char_t  ring, UShort_t sec, UShort_t strip);
     // void GetVertex(Double_t* vertexXYZ) ;
@@ -74,7 +74,8 @@ class AliFMDAnalysisTaskSharing : public AliAnalysisTask
     AliESDVertex* fEsdVertex;         // vtx info from the ESD
     Bool_t        fStatus;            // event status
     AliFMDFloatMap fLastTrackByStrip; // the last track to hit this strip
-    
+    UInt_t          fLastOrbit;
+    Float_t fEtotal;
     ClassDef(AliFMDAnalysisTaskSharing, 0); // Analysis task for FMD analysis
 };