]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
Fixes of warnings and upgrade of analysis to include Pb+Pb analysis. Background corre...
authorhdalsgaa <hdalsgaa@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 20 Mar 2009 15:10:03 +0000 (15:10 +0000)
committerhdalsgaa <hdalsgaa@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 20 Mar 2009 15:10:03 +0000 (15:10 +0000)
12 files changed:
FMD/AliFMDBackgroundCorrection.cxx
FMD/AliFMDBackgroundCorrection.h
FMD/analysis/AliFMDAnaCalibEnergyDistribution.cxx
FMD/analysis/AliFMDAnaCalibEnergyDistribution.h
FMD/analysis/AliFMDAnalysisTaskBackgroundCorrection.cxx
FMD/analysis/AliFMDAnalysisTaskCollector.cxx
FMD/analysis/AliFMDAnalysisTaskDensity.cxx
FMD/analysis/AliFMDAnalysisTaskDndeta.cxx
FMD/analysis/AliFMDAnalysisTaskSharing.cxx
FMD/analysis/RunAliEnFMDAnalysisSE.C
FMD/analysis/RunLocalFMDAnalysis.C
FMD/analysis/RunLocalFMDAnalysisSE.C

index 4ebecb4f56111da2618944d6c189a2442a590e33..0edd73a95c62cfbc15482cce66f3606deace7720 100644 (file)
@@ -71,6 +71,7 @@ AliFMDBackgroundCorrection::AliFMDInputBG::AliFMDInputBG(Bool_t hits_not_trackre
   fPrimaryArray(),
   fHitArray(),
   fHitMap(),
+  fLastTrackByStrip(),
   fPrim(0),
   fHits(0),
   fZvtxCut(0),
@@ -81,9 +82,11 @@ AliFMDBackgroundCorrection::AliFMDInputBG::AliFMDInputBG(Bool_t hits_not_trackre
   fPrevSec(-1),
   fNbinsEta(100)
 {
-  AddLoad(kTracks); 
-  if(hits_not_trackref)
+
+  if(hits_not_trackref) {
     AddLoad(kHits);
+    AddLoad(kTracks); 
+  }
   else
     AddLoad(kTrackRefs);
   AddLoad(kKinematics); 
@@ -320,12 +323,13 @@ AliFMDBackgroundCorrection::AliFMDInputBG::ProcessHit(AliFMDHit* h,
                               h->Track(),
                               h->Q());
   
+  // std::cout<<"Length: "<<h->Length()<<std::endl;
+  // std::cout<<"Is stopped "<<h->IsStop()<<"    "<<kTRUE<<std::endl;
   return retval;
 }
 //_____________________________________________________________________
 Bool_t 
-AliFMDBackgroundCorrection::AliFMDInputBG::ProcessTrackRef(AliTrackReference* tr, 
-                                                          TParticle* p) 
+AliFMDBackgroundCorrection::AliFMDInputBG::ProcessTrackRef(AliTrackReference* tr, TParticle* p) 
 {
   if(!tr)
     return kTRUE;
@@ -352,20 +356,46 @@ AliFMDBackgroundCorrection::AliFMDInputBG::ProcessEvent(UShort_t det,
   
   
   
-  if(charge !=  0 && 
+  /*  if(charge !=  0 && 
      ((nTrack != fPrevTrack) || 
-      (det != fPrevDetector) || 
-      (ring != fPrevRing)    ||
-      (sec != fPrevSec))) {
-    fHitMap.operator()(det,ring,sec,strip) = 1;
+      (det != fPrevDetector))){ || 
+      (ring != fPrevRing))){    ||
+                                 (sec != fPrevSec))) {*/
+  /*Float_t nstrips = (ring =='O' ? 256 : 512);
+  
+  
+  Float_t prevStripTrack = -1;
+  if(strip >0)
+    prevStripTrack = fLastTrackByStrip.operator()(det,ring,sec,strip-1);
+  
+  Float_t nextStripTrack = -1;
+  if(strip < (nstrips - 1))
+    nextStripTrack = fLastTrackByStrip.operator()(det,ring,sec,strip+1);
+  */
+  
+  Float_t thisStripTrack = fLastTrackByStrip.operator()(det,ring,0,0);
+  // if(nTrack == fLastTrackByStrip.operator()(det,ring,sec,strip))
+  //  std::cout<<"Track # "<<nTrack<<"  failed the cut in "<<det<<"   "<<ring<<"   "<<sec<<"  "<<strip<<std::endl;
+  // else
+  //  std::cout<<"Track # "<<nTrack<<"  passed the cut in "<<det<<"   "<<ring<<"   "<<sec<<"  "<<strip<<std::endl;
+  if(charge != 0 && nTrack != thisStripTrack) {
+    
+    fHitMap.operator()(det,ring,sec,strip) += 1;
     fHits++;
+    Float_t nstrips = (ring =='O' ? 256 : 512);
+    fLastTrackByStrip.operator()(det,ring,sec,strip) = (Float_t)nTrack;
+    if(strip >0)
+      fLastTrackByStrip.operator()(det,ring,sec,strip-1) = (Float_t)nTrack;
+    if(strip < (nstrips - 1))
+      fLastTrackByStrip.operator()(det,ring,sec,strip+1) = (Float_t)nTrack;
     
+    fPrevDetector = det;
+    fPrevRing     = ring;
+    fPrevSec      = sec;
+    fPrevTrack    = nTrack;
   }
   
-  fPrevDetector = det;
-  fPrevRing     = ring;
-  fPrevSec      = sec;
-  fPrevTrack    = nTrack;
+  
   
   return kTRUE;
 
@@ -374,6 +404,8 @@ AliFMDBackgroundCorrection::AliFMDInputBG::ProcessEvent(UShort_t det,
 //_____________________________________________________________________
 Bool_t AliFMDBackgroundCorrection::AliFMDInputBG::Init() 
 {
+  fLastTrackByStrip.Reset(-1);
+  
   fPrimaryArray.SetOwner();
   fPrimaryArray.SetName("FMD_primary");
   
@@ -457,8 +489,6 @@ Bool_t AliFMDBackgroundCorrection::AliFMDInputBG::Begin(Int_t event )
     //(charge!=0)&&(TMath::Abs(p->Vx() - vertex.At(0))<0.01)&&(TMath::Abs(p->Vy() - vertex.At(1))<0.01)&&(TMath::Abs(p->Vz() - vertex.At(2))<0.01);
     if(primary && charge !=0) {
       
-      
-      
       fPrim++;
       
       fPrimaryMapInner.Fill(eta,phi);
@@ -506,7 +536,7 @@ Bool_t AliFMDBackgroundCorrection::AliFMDInputBG::End()  {
        for(UShort_t strip = 0; strip < nstr; strip++) {
          
          if(fHitMap.operator()(det,ring,sec,strip) > 0) {
-         
+           //std::cout<<fHitMap.operator()(det,ring,sec,strip)<<std::endl;
            Double_t x,y,z;
            AliFMDGeometry* fmdgeo = AliFMDGeometry::Instance();
            fmdgeo->Detector2XYZ(det,ring,sec,strip,x,y,z);
@@ -522,7 +552,7 @@ Bool_t AliFMDBackgroundCorrection::AliFMDInputBG::End()  {
            Float_t   r     = TMath::Sqrt(TMath::Power(x,2)+TMath::Power(y,2));
            Float_t   theta = TMath::ATan2(r,z-vertex.At(2));
            Float_t   eta   = -1*TMath::Log(TMath::Tan(0.5*theta));
-           hHits->Fill(eta,phi);
+           hHits->Fill(eta,phi,fHitMap.operator()(det,ring,sec,strip));
          }
          
        }
@@ -533,6 +563,7 @@ Bool_t AliFMDBackgroundCorrection::AliFMDInputBG::End()  {
   fPrimaryMapInner.Reset();
   fPrimaryMapOuter.Reset();
   fHitMap.Reset(0);
+  fLastTrackByStrip.Reset(-1);
   
   return retval;
 }
index f51dda4d6b086d434deabc0d2bedd8fd6269445c..f30fee8e58ad676a2cfd30ee8240a5b8b8d3a79d 100644 (file)
@@ -73,6 +73,7 @@ public:
     TH2F    fPrimaryMapInner;
     TH2F    fPrimaryMapOuter;
     AliFMDFloatMap fHitMap;
+    AliFMDFloatMap fLastTrackByStrip;
     Int_t fPrim;
     Int_t fHits;
     Double_t fZvtxCut;
index 365a4fac9264b236b30a29b300487bf59187cd55..91f54dd1a415344ee228a7b376b224f4a4f63184 100644 (file)
@@ -1,43 +1,69 @@
 
 #include "AliFMDAnaCalibEnergyDistribution.h"
-
+#include "TAxis.h"
+#include <AliLog.h>
+#include <iostream>
 ClassImp(AliFMDAnaCalibEnergyDistribution)
 
 AliFMDAnaCalibEnergyDistribution::AliFMDAnaCalibEnergyDistribution() : TObject(),
-  fArray(), fIsInit(kFALSE){
-  
-  
+  fArray(), 
+  fIsInit(kFALSE),
+  fNetaBins(0),
+  fEtaMax(0),
+  fEtaMin(0){
+   
   
 }
 //____________________________________________________________________
 void AliFMDAnaCalibEnergyDistribution::Init() {
   
+  if(fNetaBins == 0)
+    AliFatal("Set Eta bins before doing Init or anything else");
+  
   fArray.SetOwner();
-  for(Int_t det = 1; det<=3;det++) {
-    TObjArray* detArray = new TObjArray();
-    fArray.AddAtAndExpand(detArray,det);
+  
+  for(Int_t i = 0; i<=fNetaBins+1; i++) {
+    TObjArray* etaArray = new TObjArray();
+    fArray.AddAtAndExpand(etaArray,i);
+    for(Int_t det = 1; det<=3;det++) {
+      TObjArray* detArray = new TObjArray();
+      etaArray->AddAtAndExpand(detArray,det);
+      
+      
+    }
+  
   }
   fIsInit = kTRUE;
 }
 
 
 //____________________________________________________________________
-TH1F* AliFMDAnaCalibEnergyDistribution::GetEnergyDistribution(Int_t det, Char_t ring) {
-
+TH1F* AliFMDAnaCalibEnergyDistribution::GetEnergyDistribution(Int_t det, Char_t ring, Float_t eta) {
+  
+  TAxis testaxis(fNetaBins,fEtaMin,fEtaMax);
+  Int_t binnumber = testaxis.FindBin(eta);
+  
   Int_t ringNumber     = (ring == 'I' ? 0 : 1);
-  TObjArray* detArray  = (TObjArray*)fArray.At(det);
+  TObjArray* etaArray  = (TObjArray*)fArray.At(binnumber); 
+  TObjArray* detArray  = (TObjArray*)etaArray->At(det); 
   TH1F* hEdist         = (TH1F*)detArray->At(ringNumber);    
   
   return hEdist;
 }
 
 //____________________________________________________________________
-void AliFMDAnaCalibEnergyDistribution::SetEnergyDistribution(Int_t det, Char_t ring, TH1F* edist) {
+void AliFMDAnaCalibEnergyDistribution::SetEnergyDistribution(Int_t det, Char_t ring, Float_t eta, TH1F* edist ) {
   
   if(!fIsInit)
     Init();
+  
+  TAxis testaxis(fNetaBins,fEtaMin,fEtaMax);
+  Int_t binnumber = testaxis.FindBin(eta);
+  
   Int_t ringNumber     = (ring == 'I' ? 0 : 1);
-  TObjArray* detArray  = (TObjArray*)fArray.At(det);
+  TObjArray* etaArray  = (TObjArray*)fArray.At(binnumber);
+  TObjArray* detArray  = (TObjArray*)etaArray->At(det);
+  
   detArray->AddAtAndExpand(edist,ringNumber);
   
 
index 3178bf00504d76ef6dcb88de82e4b24bbba019b7..10c11d8ab1a65fb8aad60cc8cdda4f97bebbd6c9 100644 (file)
@@ -11,14 +11,21 @@ class AliFMDAnaCalibEnergyDistribution : public TObject
  public:
   
   AliFMDAnaCalibEnergyDistribution();
-  void  SetEnergyDistribution(Int_t det, Char_t ring, TH1F* edist);
-  TH1F* GetEnergyDistribution(Int_t det, Char_t ring);
-  
+  void  SetNetaBins(Int_t nbins) {fNetaBins = nbins;}
+  Int_t GetNetaBins() { return fNetaBins;}
+  void  SetEtaLimits(Float_t eta_min, Float_t eta_max) {fEtaMin = eta_min; fEtaMax = eta_max;}
+  void  SetEnergyDistribution(Int_t det, Char_t ring, Float_t eta, TH1F* edist);
+  TH1F* GetEnergyDistribution(Int_t det, Char_t ring, Float_t eta);
+
   
  protected:
   void      Init();
   TObjArray fArray;
   Bool_t    fIsInit;
+  Int_t     fNetaBins;
+  Float_t   fEtaMax;
+  Float_t   fEtaMin;
+  
   ClassDef(AliFMDAnaCalibEnergyDistribution,1);
 };
 
index d62fc00a081e5ccbd802554591f220d18528e4cc..5f6946f11a33e84f0617647f7a7dd7789ec56d52 100644 (file)
@@ -177,7 +177,21 @@ void AliFMDAnalysisTaskBackgroundCorrection::Exec(Option_t */*option*/)
       TH2F* hTmp       = (TH2F*)hMultInput->Clone("hMult_from_event");
       
       hTmp->Divide(hTmp,hBg,1,1);//,"B");
-      
+      /*for(Int_t i = 1; i<=hTmp->GetNbinsX();i++) {
+       for(Int_t j = 1; j<=hTmp->GetNbinsY();j++) {
+         Float_t mult = hTmp->GetBinContent(i,j);
+         if(mult == 0) continue;
+         Float_t correction = hBg->GetBinContent(i,j);
+         
+         Float_t multcor = mult;
+         if(correction != 0)
+           multcor = multcor/correction;
+         else
+           std::cout<<"Warning! No correction for bin "<<i<<" , "<<j<<std::endl;
+         
+         hTmp->SetBinContent(i,j,multcor);
+       }
+       }*/
       hMultTotal->Add(hTmp);
       delete hTmp;
       
index d76dbc6423239321c94995cdd341474ea317aece..9e38c6f68f674545645a8f1bf5cf8c2587d5fc0f 100644 (file)
@@ -56,28 +56,33 @@ void AliFMDAnalysisTaskCollector::CreateOutputObjects()
   printf("AnalysisTaskFMD::CreateOutPutData() \n");
   
   fOutputList = new TList();//(TList*)GetOutputData(0);
+  AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
   
   fArray     = new TObjArray();
   fArray->SetName("FMD");
   fArray->SetOwner();
   TH1F* hEdist = 0;
-  for(Int_t det =1; det<=3;det++)
-    {
-      TObjArray* detArray = new TObjArray();
-      detArray->SetName(Form("FMD%d",det));
-      fArray->AddAtAndExpand(detArray,det);
-      Int_t nRings = (det==1 ? 1 : 2);
-      for(Int_t ring = 0;ring<nRings;ring++)
-       {
-         Char_t ringChar = (ring == 0 ? 'I' : 'O');
-         hEdist = new TH1F(Form("FMD%d%c",det,ringChar),Form("FMD%d%c",det,ringChar),100,0,3);
-         hEdist->SetXTitle("#Delta E / E_{MIP}");
-         fOutputList->Add(hEdist);
-         detArray->AddAtAndExpand(hEdist,ring);
-       } 
-    }
-  
-  
+  for(Int_t nEta = 0; nEta <= pars->GetNetaBins()+1; nEta++) {
+    TObjArray* etaArray = new TObjArray();
+    fArray->AddAtAndExpand(etaArray,nEta);
+    for(Int_t det =1; det<=3;det++)
+      {
+       TObjArray* detArray = new TObjArray();
+       detArray->SetName(Form("FMD%d",det));
+       etaArray->AddAtAndExpand(detArray,det);
+       Int_t nRings = (det==1 ? 1 : 2);
+       for(Int_t ring = 0;ring<nRings;ring++)
+         {
+           Char_t ringChar = (ring == 0 ? 'I' : 'O');
+           hEdist = new TH1F(Form("FMD%d%c_etabin%d",det,ringChar,nEta),Form("FMD%d%c_etabin%d",det,ringChar,nEta),200,0,6);
+           hEdist->SetXTitle("#Delta E / E_{MIP}");
+           fOutputList->Add(hEdist);
+           detArray->AddAtAndExpand(hEdist,ring);
+         } 
+      }
+    
+  }
   
   fZvtxDist  = new TH1F("ZvtxDist","Vertex distribution",100,-30,30);
   fZvtxDist->SetXTitle("z vertex");
@@ -112,8 +117,8 @@ void AliFMDAnalysisTaskCollector::Exec(Option_t */*option*/)
   Double_t vertex[3];
   
   GetVertex(vertex);
-  if(vertex[0] == 0 && vertex[1] == 0 && vertex[2] == 0)
-    return;
+  //if(vertex[0] == 0 && vertex[1] == 0 && vertex[2] == 0)
+  //  return;
   fZvtxDist->Fill(vertex[2]);
   
   if(TMath::Abs(vertex[2]) > pars->GetVtxCutZ())
@@ -122,18 +127,29 @@ void AliFMDAnalysisTaskCollector::Exec(Option_t */*option*/)
   AliESDFMD* fmd = fESD->GetFMDData();
   if (!fmd) return;
   
+  
+  
   for(UShort_t det=1;det<=3;det++) {
-    TObjArray* detArray = (TObjArray*)fArray->At(det);
+    
+    
     Int_t nRings = (det==1 ? 1 : 2);
     for (UShort_t ir = 0; ir < nRings; ir++) {
-      TH1F* Edist = (TH1F*)detArray->At(ir);
+  
       Char_t   ring = (ir == 0 ? 'I' : 'O');
       UShort_t nsec = (ir == 0 ? 20  : 40);
       UShort_t nstr = (ir == 0 ? 512 : 256);
+      TH2F* hBg = pars->GetBackgroundCorrection(det,ring,0);
+      
       for(UShort_t sec =0; sec < nsec;  sec++)  {
        for(UShort_t strip = 0; strip < nstr; strip++) {
+         Float_t eta = fmd->Eta(det,ring,sec,strip);
+         Int_t nEta = hBg->GetXaxis()->FindBin(eta);
+        
+         TObjArray* etaArray = (TObjArray*)fArray->At(nEta);
+         TObjArray* detArray = (TObjArray*)etaArray->At(det);
+         TH1F* Edist = (TH1F*)detArray->At(ir);
          Float_t mult = fmd->Multiplicity(det,ring,sec,strip);
-         if(mult == AliESDFMD::kInvalidMult) continue;
+         if(mult == AliESDFMD::kInvalidMult || mult == 0) continue;
          Edist->Fill(mult);
          
        }
index 67db948203a8aed25e0e61ef208f0e47a996d694..c10a16729795a619bdb8ea72367e824ceddba347 100644 (file)
@@ -170,20 +170,48 @@ void AliFMDAnalysisTaskDensity::Exec(Option_t */*option*/)
       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);
-         Float_t mult_cut = 0.1;
+         Float_t eta = fESD->Eta(det,ring,sec,strip);
+         Float_t mult_cut = 0.2;
          if(mult == 0 || mult == AliESDFMD::kInvalidMult) continue;
          //Particle number cut goes here...
          Float_t nParticles = 0;
-         if(mult > mult_cut)
-           nParticles = 1;
-         Float_t eta = fESD->Eta(det,ring,sec,strip);
+         if(fESD->GetUniqueID() == kTRUE) {
+           //proton + proton
+           if(mult > mult_cut) 
+             nParticles = 1; 
+         }
+         else {
+           
+           //Pb+Pb
+           Float_t mpv   = pars->GetMPV(det,ring,eta);
+           Float_t sigma = pars->GetSigma(det,ring,eta);
+           Float_t alpha = pars->Get2MIPWeight(det,ring,eta);
+           Float_t beta  = pars->Get3MIPWeight(det,ring,eta);
+           
+           Float_t sumCor = TMath::Landau(mult,mpv,sigma,kTRUE)+
+             alpha*TMath::Landau(mult,2*mpv+2*sigma*TMath::Log(2),2*sigma,kTRUE)+
+             beta*TMath::Landau(mult,3*mpv+3*sigma*TMath::Log(3),3*sigma,kTRUE);
+           Float_t weight = TMath::Landau(mult,mpv,sigma,kTRUE)+
+             2*alpha*TMath::Landau(mult,2*mpv+2*sigma*TMath::Log(2),2*sigma,kTRUE)+
+             3*beta*TMath::Landau(mult,3*mpv+3*sigma*TMath::Log(3),3*sigma,kTRUE);
+           
+           
+           if(mult > 0){//mult_cut) {
+             if(sumCor) nParticles = weight / sumCor;
+             else nParticles = 1;
+           }
+           //std::cout<<sumCor<<"    "<<weight<<"    "<<"    "<<mult<<"  "<<nParticles<<std::endl;
+           
+         }
+         
+         
          Double_t x,y,z;
          geo->Detector2XYZ(det,ring,sec,strip,x,y,z);
          Float_t phi = TMath::ATan2(y,x);
          if(phi<0)
            phi = phi+2*TMath::Pi();
          Float_t correction = GetAcceptanceCorrection(ring,strip);
-         if(correction) mult = mult / correction;
+         if(correction) nParticles = nParticles / correction;
          hMult->Fill(eta,phi,nParticles);
          
          
index acf4104928530507e974667b7c1177cbf7d77ca6..7d6491b76a8766a8973d778b06a55d3195ae983e 100644 (file)
@@ -22,7 +22,8 @@
 #include "AliFMDAnaParameters.h"
 #include "AliFMDGeometry.h"
 #include "AliGenEventHeader.h"
-
+#include "TDatabasePDG.h"
+#include "TParticlePDG.h"
 ClassImp(AliFMDAnalysisTaskDndeta)
 
 
@@ -153,7 +154,6 @@ void AliFMDAnalysisTaskDndeta::Exec(Option_t */*option*/)
 {
   Int_t vtxbin   = fVertexString->GetString().Atoi();
   fNevents.Fill(vtxbin);
-  
   for(UShort_t det=1;det<=3;det++) {
     //TObjArray* detInputArray = (TObjArray*)fInputArray->At(det);
     TObjArray* detArray = (TObjArray*)fArray.At(det);
@@ -222,11 +222,13 @@ void AliFMDAnalysisTaskDndeta::ProcessPrimary() {
       continue;
     if(TMath::Abs(particle->Zv()) > pars->GetVtxCutZ())
       continue;
+    
     if(stack->IsPhysicalPrimary(i) && particle->Charge() != 0) {
       hPrimary->Fill(particle->Eta());
       Double_t delta           = 2*pars->GetVtxCutZ()/pars->GetNvtxBins();
       Double_t vertexBinDouble = (particle->Zv() + pars->GetVtxCutZ()) / delta;
       Int_t    vertexBin       = (Int_t)vertexBinDouble;
+
       TH1F* hPrimVtxBin = (TH1F*)fOutputList->FindObject(Form("primmult_vtxbin%d",vertexBin));
       hPrimVtxBin->Fill(particle->Eta());
       if(firstTrack) {
index d96d79ca1820a3beff45bb9518558f6493c08f42..ddb384ee5c62ad17ad8687d217fc737672bbfe52 100644 (file)
@@ -15,6 +15,7 @@
 #include "AliMCEventHandler.h"
 #include "AliStack.h"
 #include "AliESDVertex.h"
+#include "AliMultiplicity.h"
 #include "AliFMDAnaParameters.h"
 #include "AliFMDParameters.h"
 
@@ -25,6 +26,8 @@ AliFMDAnalysisTaskSharing::AliFMDAnalysisTaskSharing()
 : fDebug(0),
   fESD(0x0),
   foutputESDFMD(),
+  fEnergy(0),
+  fNstrips(0),
   fSharedThis(kFALSE),
   fSharedPrev(kFALSE),
   fDiagList(),
@@ -118,11 +121,18 @@ void AliFMDAnalysisTaskSharing::Exec(Option_t */*option*/)
   fEsdVertex->SetXYZ(vertex);
   
   if(vertex[0] == 0 && vertex[1] == 0 && vertex[2] == 0) {
+  
     fStatus = kFALSE;
     return;
   }
   else
     fStatus = kTRUE;
+  const AliMultiplicity* testmult = fESD->GetMultiplicity();
+  
+  Int_t nTrackLets = testmult->GetNumberOfTracklets();
+  
+  if(nTrackLets < 1000) foutputESDFMD->SetUniqueID(kTRUE);
+  else foutputESDFMD->SetUniqueID(kFALSE);
   
   AliESDFMD* fmd = fESD->GetFMDData();
   
@@ -146,11 +156,11 @@ void AliFMDAnalysisTaskSharing::Exec(Option_t */*option*/)
          foutputESDFMD->SetMultiplicity(det,ring,sec,strip,0.);
          Float_t mult = fmd->Multiplicity(det,ring,sec,strip);
          
-         
          if(mult == AliESDFMD::kInvalidMult || mult == 0) continue;
          
          Double_t eta  = fmd->Eta(det,ring,sec,strip);//EtaFromStrip(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(fmd->Eta(det,ring,sec,strip)));
@@ -170,6 +180,7 @@ void AliFMDAnalysisTaskSharing::Exec(Option_t */*option*/)
            }
          
          Float_t merged_energy = GetMultiplicityOfStrip(mult,eta,Eprev,Enext,det,ring,sec,strip);
+
          if(merged_energy > 0 )
            nHits++;
          foutputESDFMD->SetMultiplicity(det,ring,sec,strip,merged_energy);
@@ -180,8 +191,6 @@ void AliFMDAnalysisTaskSharing::Exec(Option_t */*option*/)
     }
   }
   
-  //std::cout<<fESD->GetEventNumberInFile()<<"    "<<nHits<<"   "<<vertex[2]<<std::endl;
-  
   if(fStandalone) {
     PostData(0, foutputESDFMD); 
     PostData(1, fEsdVertex); 
@@ -196,34 +205,30 @@ Float_t AliFMDAnalysisTaskSharing::GetMultiplicityOfStrip(Float_t mult,
                                                          Float_t Enext,
                                                          UShort_t   det,
                                                          Char_t  ring,
-                                                         UShort_t sec,
-                                                         UShort_t strip) {
+                                                         UShort_t /*sec*/,
+                                                         UShort_t /*strip*/) {
   AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
-  AliFMDParameters* recopars = AliFMDParameters::Instance();
   Float_t merged_energy = 0;
   Float_t nParticles = 0;
   Float_t cutLow  = 0.2;
-  //  Float_t gain = recopars->GetPulseGain(det,ring,sec,strip);
-  //Float_t pw   = recopars->GetPedestalWidth(det,ring,sec,strip);
-  //Float_t cutLow  =  (4*pw)/(gain*recopars->GetDACPerMIP());
-  //std::cout<<cutLow<<"    "<<gain<<"    "<<recopars->GetEdepMip()<<std::endl;
-  Float_t cutHigh = pars->GetMPV(det,ring) -pars->GetSigma(det,ring);
-  Float_t cutPart = pars->GetMPV(det,ring) - 5*pars->GetSigma(det,ring);
-  // Float_t cutPart = cutLow;//0.33*pars->GetMPV(det,ring);
-  
-  // if(ring == 'O')
-  //  cutPart = pars->GetMPV(det,ring) - 2*pars->GetSigma(det,ring);
+  Float_t cutHigh = pars->GetMPV(det,ring,eta) -1*pars->GetSigma(det,ring,eta);
+  // Float_t cutPart = pars->GetMPV(det,ring,eta) - 5*pars->GetSigma(det,ring,eta);
   Float_t Etotal  = mult;
-  //std::cout<<mult<<"  "<<strip<<std::endl;
   
-  if(mult > 0 ) {
-    fEnergy = fEnergy + mult;
-    fNstrips++;
-  }
-  if((Enext <0.01 && fEnergy >0) || fNstrips >4 ) {
+  //if(mult > 5)
+  //  std::cout<<mult<<"    "<<det<<"    "<<ring<<"   "<<sec<<"    "<<strip<<std::endl;
+  
+  if(foutputESDFMD->GetUniqueID() == kTRUE) {
+    
+    if(mult > cutLow ) {
+      fEnergy = fEnergy + mult;
+      fNstrips++;
+    }
+  if((Enext <0.01 && fEnergy >0) || fNstrips >2 ) {
     
           
-    if((fEnergy*TMath::Cos(Eta2Theta(eta))) > cutPart || fNstrips > 1) {
+    //if((fEnergy*TMath::Cos(Eta2Theta(eta))) > cutPart || fNstrips > 1) {
       nParticles = 1;
       merged_energy = fEnergy*TMath::Cos(Eta2Theta(eta));
       TH1F* hEdist = (TH1F*)fDiagList.FindObject(Form("Edist_after_sharing_FMD%d%c",det,ring));
@@ -232,7 +237,7 @@ Float_t AliFMDAnalysisTaskSharing::GetMultiplicityOfStrip(Float_t mult,
       hNstrips->Fill(fNstrips);
       //  std::cout<<Form("Merged signals %f %f %f into %f , %f in strip %d, sec %d, ring %c, det %d",Eprev, mult, Enext, fEnergy/TMath::Cos(Eta2Theta(eta)),fEnergy,strip,sec,ring,det )<<std::endl;
       
-    }
+      // }
     // else
     //std::cout<<Form("NO HIT  for  %f %f %f into %f , %f in strip %d, sec %d, ring %c, det %d, cuts %f , %f",Eprev, mult, Enext, fEnergy/TMath::Cos(Eta2Theta(eta)),fEnergy,strip,sec,ring,det,cutPart,cutHigh )<<std::endl;
     
@@ -243,42 +248,22 @@ Float_t AliFMDAnalysisTaskSharing::GetMultiplicityOfStrip(Float_t mult,
   
   return 0;
   
-  
-  /*
+  }
+  else {
+     
   if(fSharedThis) {
     fSharedThis      = kFALSE;
     fSharedPrev      = kTRUE;
     return 0.;
   }
   
-  if(mult < cutLow)
-     return 0;
-     
-  
-  //if(Etotal < 0.33*pars->GetMPV(det,ring)) {
-  //  fSharedThis      = kFALSE;
-  //  fSharedPrev      = kFALSE;
-  //  return 0.; 
- // }
-  
-  //Experimental cut 
-  if(mult<Enext && Enext>cutHigh)
-    {
-      fSharedThis      = kFALSE;
-      fSharedPrev      = kFALSE;
-      return 0;
-    }
-   //Experimental cut 2
-  if(mult>Enext && mult > cutHigh)
-    {
-      fSharedThis      = kTRUE;
-      fSharedPrev      = kFALSE;
-      nParticles = 1;
-    }
-  
-  
+  if(mult < cutLow) {
+    fSharedThis      = kFALSE;
+    fSharedPrev      = kFALSE;
+    return 0;
+  }
   
-     if(Eprev > cutLow && Eprev < cutHigh && !fSharedPrev ) {
+  if(Eprev > cutLow && Eprev < cutHigh && !fSharedPrev ) {
     Etotal += Eprev;
   }
   
@@ -286,24 +271,26 @@ Float_t AliFMDAnalysisTaskSharing::GetMultiplicityOfStrip(Float_t mult,
     Etotal += Enext;
     fSharedThis      = kTRUE;
   }
-      TH1F* hEdist = (TH1F*)fDiagList.FindObject(Form("Edist_after_sharing_FMD%d%c",det,ring));
-     hEdist->Fill(Etotal);
+  TH1F* hEdist = (TH1F*)fDiagList.FindObject(Form("Edist_after_sharing_FMD%d%c",det,ring));
+  hEdist->Fill(Etotal);
   
-     Etotal = Etotal*TMath::Cos(Eta2Theta(eta));
-     if(Etotal > cutPart) {
-    nParticles = 1;
+  Etotal = Etotal*TMath::Cos(Eta2Theta(eta));
+  if(Etotal > 0) {
+    merged_energy = Etotal;
     fSharedPrev      = kTRUE;
-    // std::cout<<Form("Merged signals %f %f %f into %f , %f in strip %d, sec %d, ring %c, det %d",Eprev, mult, Enext, Etotal/TMath::Cos(Eta2Theta(eta)),Etotal,strip,sec,ring,det )<<std::endl;
+    //if(det == 3 && ring =='I')
+    //  std::cout<<Form("Merged signals %f %f %f into %f , %f in strip %d, sec %d, ring %c, det %d",Eprev, mult, Enext, Etotal/TMath::Cos(Eta2Theta(eta)),Etotal,strip,sec,ring,det )<<std::endl;
   }
-  else {
-    if(Etotal > 0)
-      //std::cout<<Form("NO HIT  for  %f %f %f into %f , %f in strip %d, sec %d, ring %c, det %d, cuts %f , %f",Eprev, mult, Enext, Etotal/TMath::Cos(Eta2Theta(eta)),Etotal,strip,sec,ring,det,cutPart,cutHigh )<<std::endl;
+    else{// if(Etotal > 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",Eprev, mult, Enext, Etotal/TMath::Cos(Eta2Theta(eta)),Etotal,strip,sec,ring,det,cutPart,cutHigh )<<std::endl;
     fSharedThis      = kFALSE;
     fSharedPrev      = kFALSE;
   }
+  // merged_energy = mult;
   
-     return nParticles
-     */
+  return merged_energy
+  }  
 }
 //_____________________________________________________________________
 void AliFMDAnalysisTaskSharing::GetVertex(Double_t* vertexXYZ) 
@@ -340,7 +327,7 @@ Float_t AliFMDAnalysisTaskSharing::Eta2Theta(Float_t eta) {
   if(eta < 0)
     theta = theta-TMath::Pi();
   
-  // std::cout<<"From eta2Theta: "<<theta<<std::endl;
+  std::cout<<"From eta2Theta: "<<theta<<"   "<<eta<<std::endl;
   return theta;
   
 
index 8aed91e13928dbcc000da49e0bab28bb0dcb6b9e..2b2ffa0fd561bb4afc7dea5008461ab6ffe2760f 100644 (file)
@@ -63,7 +63,7 @@ void RunAliEnFMDAnalysisSE(const Char_t* collectionName="collection.xml", const
   timer.Start();
   if (mgr->InitAnalysis()) {
     mgr->PrintStatus();
-    mgr->StartAnalysis("local",chain, 3000);
+    mgr->StartAnalysis("local",chain, 2000);
   }   
   timer.Stop();
   timer.Print();
@@ -83,7 +83,11 @@ TChain* CreateChainSingle(const char* xmlfile, const char *treeName="esdTree")
 
   TChain* chain = new TChain(treeName);
   myCollection->Reset() ;
-  while ( myCollection->Next() ) chain->Add(myCollection->GetTURL("")) ;
+  Int_t nFiles = 0;
+  while ( myCollection->Next() && nFiles <20) { 
+    chain->Add(myCollection->GetTURL("")) ;
+    nFiles++;
+  }
   chain->ls();
   return chain;
 }
index e6e12c884e9b36222634a2b5b4e5e4392a459d7d..4b2281e81154ae8833e6fd5d8658f2f6417ecdb9 100644 (file)
@@ -8,7 +8,7 @@ void RunLocalFMDAnalysis(const Char_t* filename= "AliESDs.root",
   
   AliCDBManager* cdb = AliCDBManager::Instance();
   cdb->SetDefaultStorage(cdbPath);
-  cdb->SetSpecificStorage("FMD/*","local://$ALICE_ROOT");
+  // cdb->SetSpecificStorage("FMD/*","local://$ALICE_ROOT");
   cdb->SetRun(0);
   //AliFMDParameters* recopars = AliFMDParameters::Instance();
   // recopars->Init();
@@ -83,7 +83,7 @@ void RunLocalFMDAnalysis(const Char_t* filename= "AliESDs.root",
   TStopwatch timer;
   timer.Start();
   cout<<"Executing analysis"<<endl;
-  mgr->StartAnalysis("local",chain,1000);
+  mgr->StartAnalysis("local",chain);
   timer.Stop();
   timer.Print();
 }
index 47fdb6547d9ca613ae638b5418a18cd3b66b6ca5..0914202ab3c29927d4f4d9b38940ca725c3a7710 100644 (file)
@@ -54,8 +54,8 @@ void RunLocalFMDAnalysisSE(const Char_t* filename= "AliESDs.root", const Char_t*
   
   AliCDBManager* cdb = AliCDBManager::Instance();
   cdb->SetDefaultStorage(cdbPath);
-  cdb->SetSpecificStorage("FMD/*","local://$ALICE_ROOT");
-  cdb->SetRun(0);
+  // cdb->SetSpecificStorage("FMD/*","local://$ALICE_ROOT");
+   cdb->SetRun(0);
   
   AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
   pars->Init();