#include "AliVertexerTracks.h"
#include "AliCDBManager.h"
#include "AliGeomManager.h"
+#include "AliGRPManager.h"
+#include "AliITSDetTypeRec.h"
#include "AliITSVertexer3D.h"
#include "AliITSVertexerZ.h"
#include "AliESDVertex.h"
+#include "AliESDEvent.h"
#include "AliRun.h"
#include "AliRunLoader.h"
#include "AliGenEventHeader.h"
/* $Id$ */
-Bool_t DoVerticesSPD(Int_t optdebug=1){
+Bool_t DoVerticesSPD(Int_t pileupalgo=1, Int_t optdebug=0){
TFile *fint = new TFile("VertexSPD.root","recreate");
- TNtuple *nt = new TNtuple("ntvert","vertices","xtrue:ytrue:ztrue:zZ:zdiffZ:zerrZ:ntrksZ:x3D:xdiff3D:xerr3D:y3D:ydiff3D:yerr3D:z3D:zdiff3D:zerr3D:ntrks3D:dndy:ntrklets:nrp1:ptyp:is3D");
+ TNtuple *nt = new TNtuple("ntvert","vertices","xtrue:ytrue:ztrue:zZ:zdiffZ:zerrZ:ntrksZ:x3D:xdiff3D:xerr3D:y3D:ydiff3D:yerr3D:z3D:zdiff3D:zerr3D:ntrks3D:dndy:ntrklets:nrp1:ptyp:is3D:isTriggered");
if (gClassTable->GetID("AliRun") < 0) {
gInterpreter->ExecuteMacro("loadlibs.C");
AliCDBManager* man = AliCDBManager::Instance();
if (!man->IsDefaultStorageSet()) {
printf("Setting a local default storage and run number 0\n");
- man->SetDefaultStorage("local://$ALICE_ROOT");
+ man->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
man->SetRun(0);
}
else {
if(!gGeoManager){
AliGeomManager::LoadGeometry("geometry.root");
}
+ AliGeomManager::ApplyAlignObjsFromCDB("ITS");
AliRunLoader* runLoader = AliRunLoader::Open("galice.root");
if (!runLoader) {
AliITSLoader* ITSloader = (AliITSLoader*) runLoader->GetLoader("ITSLoader");
ITSloader->LoadRecPoints("read");
- AliMagF * magf = gAlice->Field();
- AliTracker::SetFieldMap(magf,kTRUE);
- if(optdebug) printf("MagneticField=%f\n",AliTracker::GetBz());
-
Int_t totev=runLoader->GetNumberOfEvents();
if(optdebug) printf("Number of events= %d\n",totev);
-
+ TFile* esdFile = TFile::Open("AliESDs.root");
+ if (!esdFile || !esdFile->IsOpen()) {
+ printf("Error in opening ESD file");
+ return kFALSE;
+ }
+
+ AliESDEvent * esd = new AliESDEvent;
+ TTree* tree = (TTree*) esdFile->Get("esdTree");
+ if (!tree) {
+ printf("Error: no ESD tree found");
+ return kFALSE;
+ }
+ esd->ReadFromTree(tree);
+ AliGRPManager * grpMan=new AliGRPManager();
+ grpMan->ReadGRPEntry();
+ grpMan->SetMagField();
+ printf("Magnetic field set to %f T\n",AliTracker::GetBz()/10.);
+
+ AliITSDetTypeRec* detTypeRec = new AliITSDetTypeRec();
Double_t xnom=0.,ynom=0.;
-
AliITSVertexerZ *vertz = new AliITSVertexerZ(xnom,ynom);
vertz->Init("default");
AliITSVertexer3D *vert3d = new AliITSVertexer3D();
vert3d->Init("default");
vert3d->SetWideFiducialRegion(40.,1.);
+ vert3d->SetPileupAlgo(pileupalgo);
vert3d->PrintStatus();
+ vertz->SetDetTypeRec(detTypeRec);
+ vert3d->SetDetTypeRec(detTypeRec);
/* uncomment these lines to use diamond constrain */
// Double_t posdiam[3]={0.03,0.1,0.};
/* end lines to be uncommented to use diamond constrain */
Int_t goodz=0,good3d=0;
+ // Trigger mask
+ ULong64_t spdFO = (1 << 14);
+ ULong64_t v0left = (1 << 11);
+ ULong64_t v0right = (1 << 12);
for (Int_t iEvent = 0; iEvent < totev; iEvent++) {
TArrayF mcVertex(3);
}
if(optdebug) printf(" dNch/dy = %f\n",dNchdy);
- TTree* cltree = ITSloader->TreeR();
+ tree->GetEvent(iEvent);
+ TTree* cltree = ITSloader->TreeR();
+ ULong64_t triggerMask=esd->GetTriggerMask();
+ // MB1: SPDFO || V0L || V0R
+ Bool_t eventTriggered = (triggerMask & spdFO || ((triggerMask & v0left) || (triggerMask & v0right)));
AliESDVertex* vtxz = vertz->FindVertexForCurrentEvent(cltree);
AliESDVertex* vtx3d = vert3d->FindVertexForCurrentEvent(cltree);
AliMultiplicity *alimult = vert3d->GetMultiplicity();
ntrklets=alimult->GetNumberOfTracklets() ;
nrecp1=ntrklets+alimult->GetNumberOfSingleClusters();
}
- TString tit=vtx3d->GetTitle();
- Bool_t is3d=kFALSE;
- if(tit.Contains("3D")) is3d=kTRUE;
TDirectory *current = gDirectory;
fint->cd();
- Float_t xnt[22];
+ Float_t xnt[23];
Double_t zz = 0.;
Double_t zresz = 0.;
Double_t zdiff3d = 0.;
Int_t ntrk3d = 0;
+ Bool_t is3d=kFALSE;
if(vtx3d){
+ if(vtx3d->IsFromVertexer3D()) is3d=kTRUE;
ntrk3d = vtx3d->GetNContributors();
if(is3d && ntrk3d>0)good3d++;
x3d = vtx3d->GetXv();
xnt[19]=float(nrecp1);
xnt[20]=float(ptype);
xnt[21]=float(is3d);
+ xnt[22]=float(eventTriggered);
nt->Fill(xnt);
current->cd();