#include <TFile.h>
#include <TTree.h>
-#include <TMath.h>
#include <TList.h>
#include <TString.h>
#include <TVector3.h>
#include <TClonesArray.h>
+#include <TGeoMatrix.h>
#include "AliLog.h"
-#include "AliESD.h"
+#include "AliESDEvent.h"
#include "AliESDtrack.h"
-#include "AliKalmanTrack.h"
+#include "AliESDCaloCluster.h"
#include "AliEMCALRecPoint.h"
#include "AliRunLoader.h"
#include "AliEMCALTrack.h"
#include "AliEMCALLoader.h"
+#include "AliEMCALGeometry.h"
+#include "AliEMCALReconstructor.h"
+#include "AliEMCALRecParam.h"
+#include "AliCDBEntry.h"
+#include "AliCDBManager.h"
+#include "AliEMCALReconstructor.h"
+
#include "AliEMCALTracker.h"
ClassImp(AliEMCALTracker)
+
//
//------------------------------------------------------------------------------
//
fCutAlphaMin(-200.0),
fCutAlphaMax(200.0),
fCutAngle(100.0),
- fMaxDist(100.0),
+ fMaxDist(10.0),
+ fCutNITS(3.0),
+ fCutNTPC(20.0),
fRho(1.0),
fX0(1.0),
fTracks(0),
fClusters(0),
- fMatches(0)
+ fMatches(0),
+ fGeom(0)
{
//
// Default constructor.
// and all collections to NULL.
// Output file name is set to a default value.
//
+ InitParameters();
}
//
//------------------------------------------------------------------------------
fCutAlphaMax(copy.fCutAlphaMax),
fCutAngle(copy.fCutAngle),
fMaxDist(copy.fMaxDist),
+ fCutNITS(copy.fCutNITS),
+ fCutNTPC(copy.fCutNTPC),
fRho(copy.fRho),
fX0(copy.fX0),
fTracks((TObjArray*)copy.fTracks->Clone()),
fClusters((TObjArray*)copy.fClusters->Clone()),
- fMatches((TList*)copy.fMatches->Clone())
+ fMatches((TList*)copy.fMatches->Clone()),
+ fGeom(copy.fGeom)
{
//
// Copy constructor
fCutAlphaMax = copy.fCutAlphaMax;
fCutAngle = copy.fCutAngle;
fMaxDist = copy.fMaxDist;
+ fCutNITS = copy.fCutNITS;
+ fCutNTPC = copy.fCutNTPC;
fTracks = (TObjArray*)copy.fTracks->Clone();
fClusters = (TObjArray*)copy.fClusters->Clone();
fMatches = (TList*)copy.fMatches->Clone();
+ fGeom = copy.fGeom;
+
return (*this);
}
//
//------------------------------------------------------------------------------
//
+void AliEMCALTracker::InitParameters()
+{
+ //
+ // Retrieve initialization parameters
+ //
+
+ // Check if the instance of AliEMCALRecParam exists,
+ const AliEMCALRecParam* recParam = AliEMCALReconstructor::GetRecParam();
+
+ if(!recParam){
+ AliFatal("Reconstruction parameters for EMCAL not set!");
+ }
+ else{
+ fCutX = recParam->GetTrkCutX();
+ fCutY = recParam->GetTrkCutY();
+ fCutZ = recParam->GetTrkCutZ();
+ fMaxDist = recParam->GetTrkCutR();
+ fCutAngle = recParam->GetTrkCutAngle();
+ fCutAlphaMin = recParam->GetTrkCutAlphaMin();
+ fCutAlphaMax = recParam->GetTrkCutAlphaMax();
+ fCutNITS = recParam->GetTrkCutNITS();
+ fCutNTPC = recParam->GetTrkCutNTPC();
+ }
+
+}
+//
+//------------------------------------------------------------------------------
+//
+TTree* AliEMCALTracker::SearchTrueMatches()
+{
+ //Search through the list of
+ //track match candidates and clusters
+ //and look for true matches
+ //
+ //
+ if (!fClusters) return 0;
+ if (fClusters->IsEmpty()) return 0;
+ if (!fTracks) return 0;
+ if (fTracks->IsEmpty()) return 0;
+
+ TTree *outTree = new TTree("tree", "True matches from event");
+ Int_t indexT, indexC, label;
+ outTree->Branch("indexC", &indexC, "indexC/I");
+ outTree->Branch("indexT", &indexT, "indexT/I");
+ outTree->Branch("label", &label , "label/I");
+
+ Double_t dist=0.;
+ Int_t ic, nClusters = (Int_t)fClusters->GetEntries();
+ Int_t it, nTracks = fTracks->GetEntries();
+
+ for (ic = 0; ic < nClusters; ic++) {
+ AliEMCALMatchCluster *cluster = (AliEMCALMatchCluster*)fClusters->At(ic);
+ label = cluster->Label();
+ indexC = cluster->Index();
+ for (it = 0; it < nTracks; it++) {
+ AliEMCALTrack *track = (AliEMCALTrack*)fTracks->At(it);
+ if (TMath::Abs(track->GetSeedLabel()) != label) continue;
+ dist = CheckPair(track, cluster);
+ if (dist <= fMaxDist) {
+ indexT = track->GetSeedIndex();
+ outTree->Fill();
+ }
+ }
+ }
+
+ return outTree;
+}
+//
+//------------------------------------------------------------------------------
+//
void AliEMCALTracker::Clear(Option_t* option)
{
//
// Clearing method
- // Clears all specified arrays and the containers themselves.
+ // Deletes all objects in arrays and the arrays themselves
//
TString opt(option);
}
if (fTracks != 0x0 && clearTracks) {
- if (!fTracks->IsEmpty()) fTracks->Delete();
- delete fTracks;
- fTracks = 0;
+ fTracks->Delete();
+ delete fTracks;
+ fTracks = 0;
}
if (fClusters != 0x0 && clearClusters) {
- if (!fClusters->IsEmpty()) fClusters->Delete();
- delete fClusters;
- fClusters = 0;
+ fClusters->Delete();
+ delete fClusters;
+ fClusters = 0;
}
if (fMatches != 0x0 && clearMatches) {
- if (!fMatches->IsEmpty()) fMatches->Delete();
- delete fMatches;
- fMatches = 0;
+ fMatches->Delete();
+ delete fMatches;
+ fMatches = 0;
}
}
//
Clear("CLUSTERS");
+ cTree->SetBranchStatus("*",0); //disable all branches
+ cTree->SetBranchStatus("EMCALECARP",1); //Enable only the branch we need
+
TBranch *branch = cTree->GetBranch("EMCALECARP");
if (!branch) {
- AliError("can't get the branch with the EMCAL clusters");
+ AliError("Can't get the branch with the EMCAL clusters");
return 1;
}
- TClonesArray dummy("AliEMCALRecPoint", 10000);
- TClonesArray *clusters = &dummy;
+ TClonesArray *clusters = new TClonesArray("AliEMCALRecPoint", 1000);
branch->SetAddress(&clusters);
- Int_t nClusters = (Int_t)clusters->GetEntries();
- cTree->GetEvent(0);
- fClusters = new TObjArray(0);
+ //cTree->GetEvent(0);
+ branch->GetEntry(0);
+ Int_t nClusters = (Int_t)clusters->GetEntries();
+ if(fClusters) fClusters->Delete();
+ else fClusters = new TObjArray(0);
for (Int_t i = 0; i < nClusters; i++) {
AliEMCALRecPoint *cluster = (AliEMCALRecPoint*)clusters->At(i);
if (!cluster) continue;
- if (cluster->GetClusterType() != AliESDCaloCluster::kClusterv1) continue;
+ if (cluster->GetClusterType() != AliVCluster::kEMCALClusterv1) continue;
AliEMCALMatchCluster *matchCluster = new AliEMCALMatchCluster(i, cluster);
fClusters->AddLast(matchCluster);
}
- if (fClusters->IsEmpty()) {
- AliError("No clusters collected");
- return 1;
- }
-
- AliInfo(Form("Collected %d clusters", fClusters->GetEntries()));
+
+ branch->SetAddress(0);
+ clusters->Delete();
+ delete clusters;
+ if (fClusters->IsEmpty())
+ AliDebug(1,"No clusters collected");
+
+ AliDebug(1,Form("Collected %d clusters (RecPoints)", fClusters->GetEntries()));
return 0;
}
//
//------------------------------------------------------------------------------
//
-Int_t AliEMCALTracker::LoadClusters(AliESD *esd)
+Int_t AliEMCALTracker::LoadClusters(AliESDEvent *esd)
{
//
// Load EMCAL clusters in the form of AliESDCaloClusters,
// make sure that tracks/clusters collections are empty
Clear("CLUSTERS");
- Int_t start = esd->GetFirstEMCALCluster();
- Int_t nClustersEMC = esd->GetNumberOfEMCALClusters();
+ Int_t start = 0;
+ Int_t nClustersEMC = esd->GetNumberOfCaloClusters();
Int_t end = start + nClustersEMC;
fClusters = new TObjArray(0);
for (i = start; i < end; i++) {
AliESDCaloCluster *cluster = esd->GetCaloCluster(i);
if (!cluster) continue;
- if (cluster->GetClusterType() != AliESDCaloCluster::kClusterv1) continue;
+ if (!cluster->IsEMCAL()) continue ;
AliEMCALMatchCluster *matchCluster = new AliEMCALMatchCluster(i, cluster);
fClusters->AddLast(matchCluster);
}
- if (fClusters->IsEmpty()) {
- AliError("No clusters collected");
- return 1;
- }
+ if (fClusters->IsEmpty())
+ AliDebug(1,"No clusters collected");
- AliInfo(Form("Collected %d clusters", fClusters->GetEntries()));
+ AliDebug(1,Form("Collected %d clusters from ESD", fClusters->GetEntries()));
return 0;
}
//
//------------------------------------------------------------------------------
//
-Int_t AliEMCALTracker::LoadTracks(AliESD *esd)
+Int_t AliEMCALTracker::LoadTracks(AliESDEvent *esd)
{
//
// Load ESD tracks.
Int_t nTracks = esd->GetNumberOfTracks();
fTracks = new TObjArray(0);
- Int_t i, j;
- Bool_t isKink;
- Double_t alpha;
+ Int_t i=0, j=0;
+ Bool_t isKink=kFALSE;
+ Double_t alpha=0.;
for (i = 0; i < nTracks; i++) {
AliESDtrack *esdTrack = esd->GetTrack(i);
// set by default the value corresponding to "no match"
- esdTrack->SetEMCALcluster(-99999);
+ esdTrack->SetEMCALcluster(kUnmatched);
// if (esdTrack->GetLabel() < 0) continue;
// if (!(esdTrack->GetStatus() & AliESDtrack::kTOFout)) continue;
isKink = kFALSE;
for (j = 0; j < 3; j++) {
- if (esdTrack->GetKinkIndex(j) > 0) isKink = kTRUE;
+ if (esdTrack->GetKinkIndex(j) != 0) isKink = kTRUE;
}
if (isKink) continue;
AliEMCALTrack *track = new AliEMCALTrack(*esdTrack);
+ track->SetMass(0.13957018);
// check alpha and reject the tracks which fall outside EMCAL acceptance
alpha = track->GetAlpha() * TMath::RadToDeg();
if (alpha > -155.0 && alpha < 67.0) {
delete track;
continue;
}
+// if (!PropagateToEMCAL(track)) {
+// delete track;
+// continue;
+// }
track->SetSeedIndex(i);
track->SetSeedLabel(esdTrack->GetLabel());
fTracks->AddLast(track);
}
if (fTracks->IsEmpty()) {
- AliError("No tracks collected");
- return 1;
+ AliDebug(1,"No tracks collected");
}
- AliInfo(Form("Collected %d tracks", fTracks->GetEntries()));
+ AliDebug(1,Form("Collected %d tracks", fTracks->GetEntries()));
return 0;
}
//
//------------------------------------------------------------------------------
//
-Int_t AliEMCALTracker::PropagateBack(AliESD* esd)
+Int_t AliEMCALTracker::PropagateBack(AliESDEvent* esd)
{
//
// Main operation method.
// After executing match finding, stores in the same ESD object all infos
// and releases the object for further reconstruction steps.
//
-
+ //
+ // Note: should always return 0=OK, because otherwise all tracking
+ // is aborted for this event
+
if (!esd) {
AliError("NULL ESD passed");
return 1;
// step 1:
// if cluster array is empty, cluster are collected
// from the passed ESD, and work is done with ESDCaloClusters
- Int_t okLoadClusters;
+ Int_t okLoadClusters, nClusters;
if (!fClusters || (fClusters && fClusters->IsEmpty())) {
- AliInfo("Cluster array is empty. Loading clusters...");
okLoadClusters = LoadClusters(esd);
- if (okLoadClusters) return 2;
}
+ nClusters = fClusters->GetEntries();
// step 2:
// collect ESD tracks
- Int_t okLoadTracks = LoadTracks(esd);
+ Int_t okLoadTracks = LoadTracks(esd), nTracks;
if (okLoadTracks) return 3;
+ nTracks = fTracks->GetEntries();
// step 3:
// each track is propagated to the "R" position of each cluster.
// IF no clusters lie within the maximum allowed distance, no matches are assigned.
Int_t nMatches = CreateMatches();
if (!nMatches) {
- AliInfo("No good matches found.");
- return 4;
+ AliDebug(1,Form("#clusters = %d -- #tracks = %d --> No good matches found.", nClusters, nTracks));
+ return 0;
}
else {
- AliInfo(Form("Found %d matches", nMatches));
+ AliDebug(1,Form("#clusters = %d -- #tracks = %d --> Found %d matches.", nClusters, nTracks, nMatches));
}
// step 4:
// when more than 1 track share the same matched cluster, only the closest one is kept.
Int_t nRemoved = SolveCompetitions();
- AliInfo(Form("Removed %d duplicate matches", nRemoved));
+ AliDebug(1,Form("Removed %d duplicate matches", nRemoved));
if (nRemoved >= nMatches) {
AliError("Removed ALL matches! Check the algorithm or data. Nothing to save");
return 5;
trackID = track->GetSeedIndex();
AliESDtrack *esdTrack = esd->GetTrack(trackID);
if (!esdTrack) continue;
- if (esdTrack->GetLabel() == cluster->Label()) {
- esdTrack->SetEMCALcluster(cluster->Index());
- }
- else {
- esdTrack->SetEMCALcluster(-cluster->Index());
- }
+
+ // cut on its and tpc track hits
+ if(esdTrack->GetNcls(0)<=fCutNITS)continue;
+ if(esdTrack->GetNcls(1)<=fCutNTPC)continue;
+
+ esdTrack->SetEMCALcluster(cluster->Index());
nSaved++;
}
/*
AliESDtrack *esdTrack = esd->GetTrack(trackID);
if (!esdTrack) continue;
if (clusterID < 0) {
- esdTrack->SetEMCALcluster(-99999);
+ esdTrack->SetEMCALcluster(kUnmatched);
}
else {
AliEMCALMatchCluster *cluster = (AliEMCALMatchCluster*)fClusters->At(clusterID);
if (!cluster) continue;
- if (esdTrack->GetLabel() == cluster->Label()) {
- nGood++;
- esdTrack->SetEMCALcluster(cluster->Index());
- }
- else {
- esdTrack->SetEMCALcluster(-cluster->Index());
- }
+
+ esdTrack->SetEMCALcluster(cluster->Index());
nSaved++;
}
}
*/
- AliInfo(Form("Saved %d matches", nSaved));
+ AliDebug(1,Form("Saved %d matches", nSaved));
return 0;
}
//
// TEMP
- Bool_t isTrue = kFALSE;
+ //Bool_t isTrue = kFALSE;
// if (tr->GetSeedLabel() == cl->Label()) {
// isTrue = kTRUE;
-// cout << "TRUE MATCH!!!" << endl;
// }
// copy track into temporary variable
// check against cut on difference 'alpha - phi'
Double_t phi = TMath::ATan2(cl->Y(), cl->X());
phi = AngleDiff(phi, tr->GetAlpha());
- if (phi < fCutAlphaMin || phi > fCutAlphaMax) return distance;
+ if (phi < fCutAlphaMin || phi > fCutAlphaMax){
+ delete tr;
+ return distance;
+ }
// try to propagate to cluster radius
// (return the 'distance' value if it fails)
pos2[0] = cl->X();
pos2[1] = cl->Y();
pos2[2] = cl->Z();
- AliKalmanTrack::MeanMaterialBudget(pos1, pos2, param);
- rho = param[0];
+ MeanMaterialBudget(pos1, pos2, param);
+ rho = param[0]*param[4];
x0 = param[1];
}
else if (fTrackCorrMode == kTrackCorrFixed) {
x0 = 0.0;
}
if (fNPropSteps) {
- Int_t i;
- Double_t r;
+ Int_t i=0;
+ Double_t r=0.;
cout.setf(ios::fixed);
cout.precision(5);
- if (isTrue) cout << "Init : " << rt << ' ' << x << ' ' << y << ' ' << z << endl;
+ //if (isTrue) cout << "Init : " << rt << ' ' << x << ' ' << y << ' ' << z << endl;
for (i = 0; i < fNPropSteps; i++) {
r = rt + (rc - rt) * ((Double_t)(i+1)/(Double_t)fNPropSteps);
- if (!tr->PropagateTo(r, x0, rho)) return distance;
+ if (!tr->PropagateTo(r, x0, rho)){
+ delete tr;
+ return distance;
+ }
tr->GetXYZ(pos);
- if (isTrue) cout << "Step : " << r << ' ' << x << ' ' << y << ' ' << z << endl;
+ // if (isTrue) cout << "Step : " << r << ' ' << x << ' ' << y << ' ' << z << endl;
}
- if (isTrue) cout << "Clstr: " << rc << ' ' << cl->X() << ' ' << cl->Y() << ' ' << cl->Z() << endl;
+ //if (isTrue) cout << "Clstr: " << rc << ' ' << cl->X() << ' ' << cl->Y() << ' ' << cl->Z() << endl;
}
else {
// when no steps are used, no correction makes sense
- if (!tr->PropagateTo(rc, 0.0, 0.0)) return distance;
+ //if (!tr->PropagateTo(rc, 0.0, 0.0)) return distance;
+ if (!tr->PropagateToGlobal(cl->X(), cl->Y(), cl->Z(), 0.0, 0.0)){
+ delete tr;
+ return distance;
+ }
+ /*
+ Bool_t propOK = kFALSE;
+ cout << "START" << endl;
+ Double_t dist, rCHK, bestDist = 10000000.0;
+ for (Double_t rTMP = rc; rTMP> rc*0.95; rTMP -= 0.1) {
+ if (!tr->PropagateTo(rTMP)) continue;
+ propOK = kTRUE;
+ tr->GetXYZ(pos);
+ rCHK = TMath::Sqrt(x*x + y*y);
+ dist = TMath::Abs(rCHK - rc);
+ cout << rCHK << " vs. " << rc << endl;
+
+ if (TMath::Abs(rCHK - rc) < 0.01) break;
+ }
+ cout << "STOP" << endl;
+ if (!propOK) return distance;
+ */
}
// get global propagation of track at end of propagation
TVector3 vc(cl->X(), cl->Y(), cl->Z());
TVector3 vt(x, y, z);
Double_t angle = TMath::Abs(vc.Angle(vt)) * TMath::RadToDeg();
+ // check: where is the track?
+// Double_t r = TMath::Sqrt(pos[0]*pos[0] + pos[1]*pos[1]);
+// Double_t phiT = TMath::ATan2(pos[1], pos[0]) * TMath::RadToDeg();
+// Double_t phiC = vc.Phi() * TMath::RadToDeg();
+ //cout << "Propagated R, phiT, phiC = " << r << ' ' << phiT << ' ' << phiC << endl;
+
if (angle > fCutAngle) {
//cout << "angle" << endl;
+ delete tr;
return distance;
}
// compute differences wr to each coordinate
x -= cl->X();
- if (x > fCutX) {
+ if (TMath::Abs(x) > fCutX) {
//cout << "cut X" << endl;
+ delete tr;
return distance;
}
y -= cl->Y();
- if (y > fCutY) {
+ if (TMath::Abs(y) > fCutY) {
//cout << "cut Y" << endl;
+ delete tr;
return distance;
}
z -= cl->Z();
- if (z > fCutZ) {
+ if (TMath::Abs(z) > fCutZ) {
//cout << "cut Z" << endl;
+ delete tr;
return distance;
}
Double_t distance = 2.0 * fMaxDist;
- Double_t x0, rho;
+ Double_t x0=0., rho=0.;
if (fTrackCorrMode == kTrackCorrMMB) {
Double_t pos1[3], pos2[3], param[6];
tr->GetXYZ(pos1);
pos2[0] = cl->X();
pos2[1] = cl->Y();
pos2[2] = cl->Z();
- AliKalmanTrack::MeanMaterialBudget(pos1, pos2, param);
- rho = param[0];
+ MeanMaterialBudget(pos1, pos2, param);
+ rho = param[0]*param[4];
x0 = param[1];
}
else if (fTrackCorrMode == kTrackCorrFixed) {
TVector3 vc(cl->X(), cl->Y(), cl->Z());
// rotate the vector in order to put all clusters on a plane intersecting
// vertically the X axis; the angle depends on the sector
- Double_t clusterRot, clusterPhi = vc.Phi() * TMath::RadToDeg();
+ Double_t clusterRot=0., clusterPhi = vc.Phi() * TMath::RadToDeg();
if (clusterPhi < 0.0) clusterPhi += 360.0;
if (clusterPhi < 100.0) {
clusterRot = -90.0;
// compute the 'phi' coordinate of the intersection point to
// the EMCAL surface
Double_t x = vc.X();
- Double_t y;
+ Double_t y = 0.;
track->GetYAt(vc.X(), track->GetBz(), y);
Double_t tmp = x*TMath::Cos(track->GetAlpha()) - y*TMath::Sin(track->GetAlpha());
y = x*TMath::Sin(track->GetAlpha()) + y*TMath::Cos(track->GetAlpha());
TVector3 vdiff = vt-vc;
// compute differences wr to each coordinate
+ delete track;
if (vdiff.X() > fCutX) return distance;
if (vdiff.Y() > fCutY) return distance;
if (vdiff.Z() > fCutZ) return distance;
//
//------------------------------------------------------------------------------
//
+Double_t AliEMCALTracker::CheckPairV3
+(AliEMCALTrack *track, AliEMCALMatchCluster *cl)
+{
+ //
+ // Given a track and a cluster,
+ // propagates the first to the radius of the second.
+ // Then, checks the propagation point against all cuts.
+ // If at least a cut is not passed, a valuer equal to
+ // twice the maximum allowed distance is passed (so the value returned
+ // will not be taken into account when creating matches)
+ //
+
+ AliEMCALTrack tr(*track);
+
+ Int_t sector=-1;
+ Double_t distance = 2.0 * fMaxDist;
+ Double_t dx=0., dy=0., dz=0.;
+ Double_t phi=0., alpha=0., slope=0., tgtXnum=0., tgtXden=0., sectorWidth = 20.0 * TMath::DegToRad();
+ Double_t xcurr=0., xprop=0., param[6] = {0., 0., 0., 0., 0., 0.}, x0=0., rho=0., bz=0.;
+ Double_t x[3]= {0., 0., 0.}, x1[3]= {0., 0., 0.}, x2[3]= {0., 0., 0.};
+
+ // get initial track position
+ xcurr = tr.GetX();
+
+ // evaluate the EMCAL sector number
+ phi = cl->Phi();
+ if (phi < 0.0) phi += TMath::TwoPi();
+ sector = (Int_t)(phi / sectorWidth);
+ alpha = ((Double_t)sector + 0.5) * sectorWidth;
+ // evaluate the corresponding X for track propagation
+ slope = TMath::Tan(alpha - 0.5*TMath::Pi());
+ tgtXnum = cl->Y() - slope * cl->X();
+ tgtXden = TMath::Sqrt(1.0 + slope*slope);
+ xprop = TMath::Abs(tgtXnum / tgtXden);
+
+ // propagate by small steps
+ tr.GetXYZ(x1);
+ bz = tr.GetBz();
+ if (!tr.GetXYZAt(xprop, bz, x2)) return distance;
+ //AliKalmanTrack::MeanMaterialBudget(x1, x2, param);
+ rho = param[0]*param[4];
+ x0 = param[1];
+ if (!tr.PropagateTo(xprop, x0, rho)) return distance;
+ //if (!tr.PropagateTo(xprop, 0.0, 0.0)) return distance;
+
+ // get propagated position at the end
+ tr.GetXYZ(x);
+ dx = TMath::Abs(x[0] - cl->X());
+ dy = TMath::Abs(x[1] - cl->Y());
+ dz = TMath::Abs(x[2] - cl->Z());
+ if (dx > fCutX || dy > fCutY || dz > fCutZ) return distance;
+
+ distance = TMath::Sqrt(dx*dx + dy*dy + dz*dz);
+
+ return distance;
+}
+//
+//------------------------------------------------------------------------------
+//
+Bool_t AliEMCALTracker::PropagateToEMCAL(AliEMCALTrack *tr)
+{
+ //
+ // Propagates the track to the proximity of the EMCAL surface
+ //
+
+ Double_t xcurr=0., xtemp=0., xprop = 438.0, step = 10.0, param[6]= {0., 0., 0., 0., 0., 0.}, x0=0., rho=0., bz=0.;
+ Double_t x1[3]= {0., 0., 0.}, x2[3]= {0., 0., 0.};
+
+ // get initial track position
+ xcurr = tr->GetX();
+
+ // propagate by small steps
+ for (xtemp = xcurr + step; xtemp < xprop; xtemp += step) {
+ // to compute material budget, take current position and
+ // propagated hypothesis without energy loss
+ tr->GetXYZ(x1);
+ bz = tr->GetBz();
+ if (!tr->GetXYZAt(xtemp, bz, x2)) return kFALSE;
+ MeanMaterialBudget(x1, x2, param);
+ rho = param[0]*param[4];
+ x0 = param[1];
+ if (!tr->PropagateTo(xtemp, x0, rho)) return kFALSE;
+ }
+
+ return kTRUE;
+}
+//
+//------------------------------------------------------------------------------
+//
Int_t AliEMCALTracker::CreateMatches()
{
//
// initialize counters and indexes
Int_t count = 0;
- Int_t ic, nClusters = (Int_t)fClusters->GetEntries();
- Int_t it, nTracks = fTracks->GetEntries();
+ Int_t ic=0, nClusters = (Int_t)fClusters->GetEntries();
+ Int_t it=0, nTracks = fTracks->GetEntries();
// external loop on clusters, internal loop on tracks
- Double_t dist;
+ Double_t dist=0.;
for (ic = 0; ic < nClusters; ic++) {
- cout << "\rMatching cluster " << ic+1 << " of " << nClusters << flush;
AliEMCALMatchCluster *cluster = (AliEMCALMatchCluster*)fClusters->At(ic);
for (it = 0; it < nTracks; it++) {
AliEMCALTrack *track = (AliEMCALTrack*)fTracks->At(it);
dist = CheckPair(track, cluster);
+ //cout << dist << endl;
if (dist <= fMaxDist) {
AliEMCALMatch *candidate = new AliEMCALMatch;
candidate->SetIndexT(it);
candidate->SetIndexC(ic);
candidate->SetDistance(dist);
+ candidate->SetPt(track->GetSignedPt());
fMatches->Add(candidate);
count++;
}
}
}
- cout << endl;
-
- /*
- // loop on clusters and tracks
- Int_t icBest;
- Double_t dist, distBest;
- for (it = 0; it < nTracks; it++) {
- AliEMCALTrack *track = (AliEMCALTrack*)fTracks->At(it);
- if (!track) continue;
- icBest = -1;
- distBest = fMaxDist;
- for (ic = 0; ic < nClusters; ic++) {
- AliEMCALMatchCluster *cluster = (AliEMCALMatchCluster*)fClusters->At(ic);
- if (!cluster) continue;
- dist = CheckPair(track, cluster);
- if (dist < distBest) {
- distBest = dist;
- icBest = ic;
- }
- }
- if (icBest >= 0) {
- track->SetMatchedClusterIndex(icBest);
- track->SetMatchedClusterDist(distBest);
- count++;
- }
- else {
- track->SetMatchedClusterIndex(-1);
- }
- }
- */
return count;
}
Int_t count = 0;
// initialize flags to check repetitions
- Int_t ic, nClusters = (Int_t)fClusters->GetEntries();
- Int_t it, nTracks = fTracks->GetEntries();
+ Int_t ic=0, nClusters = (Int_t)fClusters->GetEntries();
+ Int_t it=0, nTracks = fTracks->GetEntries();
Bool_t *usedC = new Bool_t[nClusters];
Bool_t *usedT = new Bool_t[nTracks];
for (ic = 0; ic < nClusters; ic++) usedC[ic] = kFALSE;
}
}
- /*
- Int_t it1, it2, nTracks = (Int_t)fTracks->GetEntries();
- AliEMCALTrack *track1 = 0, *track2 = 0;
- for (it1 = 0; it1 < nTracks; it1++) {
- track1 = (AliEMCALTrack*)fTracks->At(it1);
- if (!track1) continue;
- if (track1->GetMatchedClusterIndex() < 0) continue;
- for (it2 = it1+1; it2 < nTracks; it2++) {
- track2 = (AliEMCALTrack*)fTracks->At(it2);
- if (!track2) continue;
- if (track2->GetMatchedClusterIndex() < 0) continue;
- if (track1->GetMatchedClusterIndex() != track2->GetMatchedClusterIndex()) continue;
- count++;
- if (track1->GetMatchedClusterDist() < track2->GetMatchedClusterDist()) {
- track2->SetMatchedClusterIndex(-1);
- }
- else if (track2->GetMatchedClusterDist() < track1->GetMatchedClusterDist()) {
- track1->SetMatchedClusterIndex(-1);
- }
- }
- }
- */
-
+ delete [] usedC;
+ delete [] usedT;
+
return count;
}
//
void AliEMCALTracker::UnloadClusters()
{
//
- // Free memory from clusters
+ // Free memory from all arrays
+ // This method is called after the local tracking step
+ // so we can safely delete everything
//
- Clear("CLUSTERS");
+ Clear();
}
+//
+//------------------------------------------------------------------------------
+//
+TVector3 AliEMCALTracker::FindExtrapolationPoint(Double_t x,Double_t y,Double_t z, AliESDtrack *track)
+{
+ //Method to determine extrapolation point of track at location x,y,z
+ AliEMCALTrack *tr = new AliEMCALTrack(*track);
+ TVector3 error(-100.,-100.,-100.);
+ if (!tr->PropagateToGlobal(x,y,z, 0.0, 0.0)) {
+ return error;
+ }
+ Double_t pos[3]= {0., 0., 0.};
+ tr->GetXYZ(pos);
+ TVector3 ExTrPos(pos[0],pos[1],pos[2]);
+ return ExTrPos;
+}
+
//
//------------------------------------------------------------------------------
//
AliEMCALTracker::AliEMCALMatchCluster::AliEMCALMatchCluster(Int_t index, AliEMCALRecPoint *recPoint)
: fIndex(index),
- fLabel(recPoint->GetPrimaryIndex()),
+ fLabel(recPoint->GetPrimaryIndex()), //wrong! fixed below
fX(0.),
fY(0.),
fZ(0.)
fX = clpos.X();
fY = clpos.Y();
fZ = clpos.Z();
+
+ //AliEMCALRecPoint stores the track labels in the parents
+ //list, sorted according to the fractional contribution to the
+ //RecPoint. The zeroth parent gave the highest contribution
+ Int_t multparent = 0;
+ Int_t *parents = recPoint->GetParents(multparent);
+ if(multparent > 0)
+ fLabel = parents[0];
+ else
+ fLabel = -1;
+
}
//
//------------------------------------------------------------------------------
//
AliEMCALTracker::AliEMCALMatchCluster::AliEMCALMatchCluster(Int_t index, AliESDCaloCluster *caloCluster)
: fIndex(index),
- fLabel(caloCluster->GetPrimaryIndex()),
+ fLabel(caloCluster->GetLabel()),
fX(0.),
fY(0.),
fZ(0.)
// Translates an AliESDCaloCluster object into the internal format.
// Index of passed cluster in its native array must be specified.
//
- Float_t clpos[3];
- caloCluster->GetGlobalPosition(clpos);
+ Float_t clpos[3]= {0., 0., 0.};
+ caloCluster->GetPosition(clpos);
fX = (Double_t)clpos[0];
fY = (Double_t)clpos[1];
AliEMCALTracker::AliEMCALMatch *that = (AliEMCALTracker::AliEMCALMatch*)obj;
- Double_t thisDist = fDistance;
- Double_t thatDist = that->GetDistance();
+ Double_t thisDist = fPt;//fDistance;
+ Double_t thatDist = that->fPt;//that->GetDistance();
if (thisDist > thatDist) return 1;
else if (thisDist < thatDist) return -1;
fCanBeSaved(kFALSE),
fIndexC(0),
fIndexT(0),
- fDistance(0.)
+ fDistance(0.),
+ fPt(0.)
{
//default constructor
fCanBeSaved(copy.fCanBeSaved),
fIndexC(copy.fIndexC),
fIndexT(copy.fIndexT),
- fDistance(copy.fDistance)
+ fDistance(copy.fDistance),
+ fPt(copy.fPt)
{
//copy ctor
}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff