// Construct histograms and displays them.
// Use the macro EditorBar.C for best access to the functionnalities
//
-//*-- Author: Y. Schutz (SUBATECH)
+//*-- Author: Y. Schutz (SUBATECH) & Gines Martinez (SUBATECH)
//////////////////////////////////////////////////////////////////////////////
// --- ROOT system ---
if (ievent==0) cout << "AnalyzeManyEvents > " << "Starting Analyzing " << endl ;
//========== Create the Clusterizer
fClu = new AliPHOSClusterizerv1() ;
- fClu->SetEmcEnergyThreshold(0.025) ;
+ fClu->SetEmcEnergyThreshold(0.05) ;
fClu->SetEmcClusteringThreshold(0.50) ;
fClu->SetPpsdEnergyThreshold (0.0000002) ;
fClu->SetPpsdClusteringThreshold(0.0000001) ;
//fhElectronPositionY->Fill(recpart. ) ;
cout << "ELECTRON" << endl;
break;
- case kNEUTRALHADRON:
+ case kNEUTRAL_HA:
fhNeutralHadronEnergy->Fill(recparticle->Energy() ) ;
//fhNeutralHadronPositionX->Fill(recpart. ) ;
//fhNeutralHadronPositionY->Fill(recpart. ) ;
cout << "NEUTRAl HADRON" << endl;
break ;
- case kNEUTRALEM:
+ case kNEUTRAL_EM:
fhNeutralEMEnergy->Fill(recparticle->Energy() ) ;
//fhNeutralEMPositionX->Fill(recpart. ) ;
//fhNeutralEMPositionY->Fill(recpart. ) ;
//cout << "NEUTRAL EM" << endl;
break ;
- case kCHARGEDHADRON:
+ case kCHARGED_HA:
fhChargedHadronEnergy->Fill(recparticle->Energy() ) ;
//fhChargedHadronPositionX->Fill(recpart. ) ;
//fhChargedHadronPositionY->Fill(recpart. ) ;
cout << "CHARGED HADRON" << endl;
break ;
- case kGAMMAHADRON:
+ case kGAMMA_HA:
fhPhotonHadronEnergy->Fill(recparticle->Energy() ) ;
//fhPhotonHadronPositionX->Fill(recpart. ) ;
//fhPhotonHadronPositionY->Fill(recpart. ) ;
fClu = new AliPHOSClusterizerv1() ;
fClu->SetEmcEnergyThreshold(0.030) ;
- fClu->SetEmcClusteringThreshold(1.0) ;
+ fClu->SetEmcClusteringThreshold(0.50) ;
fClu->SetPpsdEnergyThreshold (0.0000002) ;
fClu->SetPpsdClusteringThreshold(0.0000001) ;
fClu->SetLocalMaxCut(0.03) ;
fPID = new AliPHOSPIDv1() ;
cout << "AnalyzeOneEvent > using particle identifier " << fPID->GetName() << endl ;
-
+ //fPID->SetShowerProfileCuts(Float_t l1m, Float_t l1M, Float_t l2m, Float_t l2M) ;
+ fPID->SetShowerProfileCuts(0.7, 2.0 , 0.6 , 1.5) ;
+
//========== Creates the Reconstructioner
- fRec = new AliPHOSReconstructioner(fClu, fTrs, fPID) ;
+ fRec = new AliPHOSReconstructioner(fClu, fTrs, fPID) ;
+ fRec -> SetDebugReconstruction(kTRUE);
//=========== Connect the various Tree's for evt
case kELECTRON:
name = "ELECTRON" ;
break ;
- case kNEUTRAL:
- name = "NEUTRAL" ;
+ case kCHARGED_HA:
+ name = "CHARGED_HA" ;
break ;
- case kCHARGEDHADRON:
- name = "CHARGED HADRON" ;
+ case kNEUTRAL_HA:
+ name = "NEUTRAL_HA" ;
break ;
- case kNEUTRALHADRON:
- name = "NEUTRAL HADRON" ;
+ case kNEUTRAL_EM:
+ name = "NEUTRAL_EM" ;
break ;
- case kNEUTRALEM:
- name = "NEUTRAL EM" ;
- break ;
- case kGAMMAHADRON:
- name = "PHOTON HADRON" ;
+ case kGAMMA_HA:
+ name = "PHOTON_HA" ;
break ;
}
typedef TClonesArray FastRecParticlesList ;
-const static Int_t kUNDEFINED = -1;
-const static Int_t kGAMMA = 0 ;
-const static Int_t kELECTRON = 1 ;
-const static Int_t kNEUTRAL = 2 ;
-const static Int_t kCHARGED = 3 ;
-const static Int_t kCHARGEDHADRON = 4 ;
-const static Int_t kNEUTRALHADRON = 5 ;
-const static Int_t kNEUTRALEM = 6 ;
-const static Int_t kGAMMAHADRON = 7 ;
+const static Int_t kUNDEFINED =-1 ;
+const static Int_t kNEUTRAL_EM = 0 ;
+const static Int_t kNEUTRAL_HA = 1 ;
+const static Int_t kGAMMA = 2 ;
+const static Int_t kGAMMA_HA = 3 ;
+const static Int_t kABSURD_EM = 4 ;
+const static Int_t kABSURD_HA = 5 ;
+const static Int_t kELECTRON = 6 ;
+const static Int_t kCHARGED_HA = 7 ;
class AliPHOSFastRecParticle : public TParticle {
// Implementation version v1 of the PHOS particle identifier
// Identification is based on information from PPSD and EMC
//
-//*-- Author: Yves Schutz (SUBATECH)
+//*-- Author: Yves Schutz (SUBATECH) & Gines Martinez (SUBATECH)
// --- ROOT system ---
Int_t index = 0 ;
AliPHOSRecParticle * rp ;
Int_t type ;
+ Int_t shower_profile; // 0 narrow and 1 wide
+ Int_t cpv_detector; // 1 hit and 0 no hit
+ Int_t pc_detector; // 1 hit and 0 no hit
while ( (tracksegment = (AliPHOSTrackSegment *)next()) ) {
new( (*rpl)[index] ) AliPHOSRecParticle(tracksegment) ;
rp = (AliPHOSRecParticle *)rpl->At(index) ;
-
- // try to figure out the type of particle:
- // 1. just looking at the PPSD information
-
- if( tracksegment->GetPpsdUpRecPoint() == 0 ) { // Neutral
-
- if( tracksegment->GetPpsdLowRecPoint() == 0 ) // Neutral
- type = kNEUTRAL ;
- else { // check the shower profile
- AliPHOSEmcRecPoint * recp = tracksegment->GetEmcRecPoint() ;
- Float_t * lambda = new Float_t[2];
- recp->GetElipsAxis(lambda) ;
- if ( ( lambda[0] > fLambda1m && lambda[0] < fLambda1M ) && // shower profile cut
- ( lambda[1] > fLambda2m && lambda[1] < fLambda2M ) )
- type = kGAMMA ; // a well identified photon
- else
- type = kGAMMAHADRON ; // looks like a photon but is a hadron (most likely)
- }
- } // Neutral
- else // Charged
- type = kCHARGED ;
-
- // 2. from the shower profile analysis
- if ( type == kNEUTRAL ) {
- AliPHOSEmcRecPoint * recp = tracksegment->GetEmcRecPoint() ;
- Float_t * lambda = new Float_t[2];
- recp->GetElipsAxis(lambda) ;
- if ( ( lambda[0] > fLambda1m && lambda[0] < fLambda1M ) && // shower profile cut
- ( lambda[1] > fLambda2m && lambda[1] < fLambda2M ) )
- type = kNEUTRALEM ;
- else
- type = kNEUTRALHADRON ;
- delete lambda ;
- }
-
- // 3. from the shower dispersion
- if (type == kCHARGED) {
-
- if( tracksegment->GetEmcRecPoint()->GetDispersion() > fCutOnDispersion) // shower dispersion cut
- type = kCHARGEDHADRON ;
- // else
- // type = kELECTRON ;
- }
+ AliPHOSEmcRecPoint * recp = tracksegment->GetEmcRecPoint() ;
+ Float_t * lambda = new Float_t[2];
+ recp->GetElipsAxis(lambda) ;
+
+ // Looking at the lateral development of the shower
+ if ( ( lambda[0] > fLambda1m && lambda[0] < fLambda1M ) && // shower profile cut
+ ( lambda[1] > fLambda2m && lambda[1] < fLambda2M ) ) shower_profile = 0 ; // NARROW PROFILE
+ else shower_profile = 1 ;// WIDE PROFILE
+
+ // Looking at the photon conversion detector
+ if( tracksegment->GetPpsdLowRecPoint() == 0 ) pc_detector = 0; // No hit
+ else pc_detector = 1; // hit
+
+ // Looking at the photon conversion detector
+ if( tracksegment->GetPpsdUpRecPoint() == 0 ) cpv_detector = 0; // No hit
+ else cpv_detector = 1; // Hit
+
+ type = shower_profile + 2*pc_detector + 4*cpv_detector;
rp->SetType(type) ;
index++ ;
}
fPx = kenergy * momdir.X() ;
fPy = kenergy * momdir.Y() ;
fPz = kenergy * momdir.Z() ;
- fType = kUNDEFINED ;
+ fType = kUNDEFINED;
fE = kenergy ; // !!! all particles have mass = 0
}
// Launches the Reconstruction process in the sequence: Make the reconstructed poins (clusterize)
// Make the track segments
// Make the reconstructed particles
-
Int_t index ;
// Digit Debuging
-
-
- if (fDebugReconstruction)
- {
- cout << ">>>>>>>>>>>>>>>>>>>>>> DebugReconstruction <<<<<<<<<<<<<<<<<<<<<<<<<<" << endl ;
- cout << "DebugReconstruction>>> Digit list entries is " << dl->GetEntries() << endl ;
- AliPHOSDigit * digit;
- Bool_t calorimeter ;
- Float_t factor;
- cout << "DebugReconstruction>>> Vol Id " <<
- " Ene(MeV, KeV) " <<
- " Index " <<
- " Nprim " <<
- " Prim1 " <<
- " Prim2 " <<
- " Prim3 " << endl;
-
- for (index = 0 ; index < dl->GetEntries() ; index++) {
- digit = (AliPHOSDigit * ) dl->At(index) ;
- calorimeter = fClusterizer->IsInEmc(digit);
- if (calorimeter) factor =1000. ; else factor=1000000.;
- cout << "DebugReconstruction>>> " <<
- setw(8) << digit->GetId() << " " <<
- setw(3) << calorimeter <<
- setw(10) << factor*fClusterizer->Calibrate(digit->GetAmp()) << " " <<
- setw(6) << digit->GetIndexInList() << " " <<
- setw(5) << digit->GetNprimary() <<" " <<
- setw(5) << digit->GetPrimary(1) <<" " <<
- setw(5) << digit->GetPrimary(2) <<" " <<
- setw(5) << digit->GetPrimary(3) << endl;
- }
-
+ if (fDebugReconstruction) {
+ cout << ">>>>>>>>>>>>>>>>>>>>>> DebugReconstruction <<<<<<<<<<<<<<<<<<<<<<<<<<" << endl ;
+ cout << "DebugReconstruction>>> Digit list entries is " << dl->GetEntries() << endl ;
+ AliPHOSDigit * digit;
+ Bool_t calorimeter ;
+ Float_t factor;
+ cout << "DebugReconstruction>>> Vol Id " <<
+ " Ene(MeV, KeV) " <<
+ " Index " <<
+ " Nprim " <<
+ " Prim1 " <<
+ " Prim2 " <<
+ " Prim3 " << endl;
+ for (index = 0 ; index < dl->GetEntries() ; index++) {
+ digit = (AliPHOSDigit * ) dl->At(index) ;
+ calorimeter = fClusterizer->IsInEmc(digit);
+ if (calorimeter) factor =1000. ; else factor=1000000.;
+ cout << "DebugReconstruction>>> " <<
+ setw(8) << digit->GetId() << " " <<
+ setw(3) << calorimeter <<
+ setw(10) << factor*fClusterizer->Calibrate(digit->GetAmp()) << " " <<
+ setw(6) << digit->GetIndexInList() << " " <<
+ setw(5) << digit->GetNprimary() <<" " <<
+ setw(5) << digit->GetPrimary(1) <<" " <<
+ setw(5) << digit->GetPrimary(2) <<" " <<
+ setw(5) << digit->GetPrimary(3) << endl;
}
+
+ }
ppsdrp = (AliPHOSPpsdRecPoint * )ppsdl->At(index) ;
ppsdrp->SetIndexInList(index) ;
}
-
+
if (fDebugReconstruction)
{
cout << "DebugReconstruction>>> Cluster emc list entries is " << emccl->GetEntries() << endl ;
AliPHOSEmcRecPoint * recpoint;
cout << "DebugReconstruction>>> Module " <<
- "Ene(MeV) " <<
- "Index " <<
- "Multi " <<
- " X " <<
- " Y " <<
- " Z " <<
- " Lambda 1 " <<
- " Lambda 2 " <<
- "MaxEnergy(MeV) " <<
+ "Ene(MeV) " <<
+ "Index " <<
+ "Multi " <<
+ " X " <<
+ " Y " <<
+ " Z " <<
+ " Lambda 1 " <<
+ " Lambda 2 " <<
+ "MaxEnergy(MeV) " <<
"Nprim " <<
"Prim1 " <<
"Prim2 " <<
setw(4) << primaries[1] << " " <<
setw(4) << primaries[2] << " " << endl;
}
-
+
cout << "DebugReconstruction>>> Cluster ppsd list entries is " << ppsdl->GetEntries() << endl ;
AliPHOSPpsdRecPoint * ppsdrecpoint;
Text_t detector[4];
cout << "DebugReconstruction>>> Module " <<
- "Det " <<
- "Ene(KeV) " <<
- "Index " <<
- "Multi " <<
- " X " <<
- " Y " <<
- " Z " <<
+ "Det " <<
+ "Ene(KeV) " <<
+ "Index " <<
+ "Multi " <<
+ " X " <<
+ " Y " <<
+ " Z " <<
"Nprim " <<
"Prim1 " <<
"Prim2 " <<
setw(4) << primaries[2] << " " << endl;
}
}
-
-
+
+
if (fDebugReconstruction) cout << "DebugReconstruction>>>> Start making track segments(unfolding+tracksegments)" << endl;
fTrackSegmentMaker->MakeTrackSegments(dl, emccl, ppsdl, trsl) ;
-
+
// mark the position of the TrackSegments in the array
AliPHOSTrackSegment * trs ;
for (index = 0 ; index < trsl->GetEntries() ; index++) {
"Prim1 " <<
"Prim2 " <<
"Prim3 " << endl;
-
+
for (index = 0 ; index < trsl->GetEntries() ; index++) {
trs = (AliPHOSTrackSegment * )trsl->At(index) ;
TVector3 locpos; trs->GetPosition(locpos);
setw(4) << primaries[0] << " " <<
setw(4) << primaries[1] << " " <<
setw(4) << primaries[2] << " " << endl;
- }
-
+ }
+
}
- if (fDebugReconstruction) cout << "DebugReconstruction>>>> Start making particles" << endl;
+ if (fDebugReconstruction) cout << "DebugReconstruction>>>> Start making reconstructed particles" << endl;
fPID->MakeParticles(trsl, rpl) ;
rp = (AliPHOSRecParticle * )rpl->At(index) ;
rp->SetIndexInList(index) ;
}
+ //Debugger of RecParticles
+ if (fDebugReconstruction){
+ cout << "DebugReconstruction>>> Reconstructed particle list entries is " << rpl->GetEntries() << endl ;
+ cout << "DebugReconstruction>>> Module " <<
+ " PARTICLE " <<
+ "Ene(KeV) " <<
+ "Index " <<
+ " X " <<
+ " Y " <<
+ " Z " <<
+ "Nprim " <<
+ "Prim1 " <<
+ "Prim2 " <<
+ "Prim3 " << endl;
+ for (index = 0 ; index < rpl->GetEntries() ; index++) {
+ rp = (AliPHOSRecParticle * ) rpl->At(index) ;
+ TVector3 locpos; (rp->GetPHOSTrackSegment())->GetPosition(locpos);
+ Int_t * primaries;
+ Int_t nprimaries;
+ Text_t particle[11];
+ primaries = (rp->GetPHOSTrackSegment())->GetPrimariesEmc(nprimaries);
+ switch(rp->GetType())
+ {
+ case kNEUTRAL_EM:
+ particle = "NEUTRAL_EM";
+ break;
+ case kNEUTRAL_HA:
+ particle = "NEUTRAL_HA";
+ break;
+ case kGAMMA:
+ particle = "GAMMA ";
+ break ;
+ case kGAMMA_HA:
+ particle = "GAMMA_HA ";
+ break ;
+ case kABSURD_EM:
+ particle = "ABSURD_EM " ;
+ break ;
+ case kABSURD_HA:
+ particle = "ABSURD_HA " ;
+ break ;
+ case kELECTRON:
+ particle = "ELECTRON " ;
+ break ;
+ case kCHARGED_HA:
+ particle = "CHARGED_HA" ;
+ break ;
+ }
+
+ cout << "DebugReconstruction>>> " <<
+ setw(4) << (rp->GetPHOSTrackSegment())->GetPHOSMod() << " " <<
+ setw(15) << particle << " " <<
+ setw(9) << 1000.*(rp->GetPHOSTrackSegment())->GetEnergy() << " " <<
+ setw(3) << rp->GetIndexInList() << " " <<
+ setw(9) << locpos.X() <<" " <<
+ setw(9) << locpos.Y() <<" " <<
+ setw(9) << locpos.Z() << " " <<
+ setw(4) << nprimaries << " " <<
+ setw(4) << primaries[0] << " " <<
+ setw(4) << primaries[1] << " " <<
+ setw(4) << primaries[2] << " " << endl;
+ }
+
+ }
+
+
}
else {
ran = fRan.Rndm() ;
if( ran <= 0.9498 )
- rv = kNEUTRALEM ;
+ rv = kNEUTRAL_EM ;
else
- rv = kNEUTRALHADRON ;
+ rv = kNEUTRAL_HA ;
}
break ;
case 2112: // it's a neutron
ran = fRan.Rndm() ;
if ( ran <= 0.9998 )
- rv = kNEUTRALHADRON ;
+ rv = kNEUTRAL_HA ;
else
- rv = kNEUTRALEM ;
+ rv = kNEUTRAL_EM ;
break ;
case -2112: // it's a anti-neutron
ran = fRan.Rndm() ;
if ( ran <= 0.9984 )
- rv = kNEUTRALHADRON ;
+ rv = kNEUTRAL_HA ;
else
- rv = kNEUTRALEM ;
+ rv = kNEUTRAL_EM ;
break ;
case 11: // it's a electron
if ( ran <= 0.9996 )
rv = kELECTRON ;
else
- rv = kCHARGEDHADRON ;
+ rv = kCHARGED_HA ;
break;
case -11: // it's a positon
if ( ran <= 0.9996 )
rv = kELECTRON ;
else
- rv = kCHARGEDHADRON ;
+ rv = kCHARGED_HA ;
break;
case -1: // it's a charged
ran = fRan.Rndm() ;
if ( ran <= 0.9996 )
- rv = kCHARGEDHADRON ;
+ rv = kCHARGED_HA ;
else
rv = kGAMMA ;