* provided "as is" without express or implied warranty. *
**************************************************************************/
+/* $Id$ */
+
////////////////////////////////////
//
// MUON track reconstructor in ALICE (class renamed from AliMUONEventReconstructor)
#include <Riostream.h>
#include <TDirectory.h>
#include <TFile.h>
-#include <TMatrixD.h> //AZ
+#include <TMatrixD.h>
#include "AliMUONTrackReconstructor.h"
#include "AliMUON.h"
#include "AliMUONHitForRec.h"
#include "AliMUONTriggerTrack.h"
#include "AliMUONTriggerCircuit.h"
+#include "AliMUONTriggerCircuitNew.h"
#include "AliMUONRawCluster.h"
#include "AliMUONLocalTrigger.h"
#include "AliMUONGlobalTrigger.h"
#include "AliRun.h" // for gAlice
#include "AliRunLoader.h"
#include "AliLoader.h"
-#include "AliMUONTrackK.h" //AZ
+#include "AliMUONTrackK.h"
#include "AliMC.h"
#include "AliLog.h"
#include "AliTrackReference.h"
//************* Defaults parameters for reconstruction
const Double_t AliMUONTrackReconstructor::fgkDefaultMinBendingMomentum = 3.0;
-const Double_t AliMUONTrackReconstructor::fgkDefaultMaxBendingMomentum = 2000.0;
+const Double_t AliMUONTrackReconstructor::fgkDefaultMaxBendingMomentum = 3000.0;
const Double_t AliMUONTrackReconstructor::fgkDefaultMaxChi2 = 100.0;
const Double_t AliMUONTrackReconstructor::fgkDefaultMaxSigma2Distance = 16.0;
const Double_t AliMUONTrackReconstructor::fgkDefaultBendingResolution = 0.01;
ClassImp(AliMUONTrackReconstructor) // Class implementation in ROOT context
- //__________________________________________________________________________
-AliMUONTrackReconstructor::AliMUONTrackReconstructor(AliLoader* loader)
+//__________________________________________________________________________
+AliMUONTrackReconstructor::AliMUONTrackReconstructor(AliLoader* loader, AliMUONData* data)
: TObject()
{
// Constructor for class AliMUONTrackReconstructor
fNRecTracks = 0; // really needed or GetEntriesFast sufficient ????
// Memory allocation for the TClonesArray of hits on reconstructed tracks
// Is 100 the right size ????
- //AZ fRecTrackHitsPtr = new TClonesArray("AliMUONTrack", 100);
fRecTrackHitsPtr = new TClonesArray("AliMUONTrackHit", 100);
fNRecTrackHits = 0; // really needed or GetEntriesFast sufficient ????
fLoader = loader;
// initialize container
- fMUONData = new AliMUONData(fLoader,"MUON","MUON");
+ // fMUONData = new AliMUONData(fLoader,"MUON","MUON");
+ fMUONData = data;
return;
}
{
// To reconstruct one event
AliDebug(1,"Enter EventReconstruct");
+ ResetTracks(); //AZ
+ ResetTrackHits(); //AZ
+ ResetSegments(); //AZ
+ ResetHitsForRec(); //AZ
MakeEventToBeReconstructed();
MakeSegments();
MakeTracks();
// To reset the array and the number of HitsForRec,
// and also the number of HitsForRec
// and the index of the first HitForRec per chamber
- if (fHitsForRecPtr) fHitsForRecPtr->Clear();
+ if (fHitsForRecPtr) fHitsForRecPtr->Delete();
fNHitsForRec = 0;
for (Int_t ch = 0; ch < AliMUONConstants::NTrackingCh(); ch++)
fNHitsForRecPerChamber[ch] = fIndexOfFirstHitForRecPerChamber[ch] = 0;
AliRunLoader *runLoader = fLoader->GetRunLoader();
AliDebug(1,"Enter MakeEventToBeReconstructed");
- ResetHitsForRec();
+ //AZ ResetHitsForRec();
if (fRecTrackRefHits == 1) {
// Reconstruction from track ref. hits
// Back to the signal file
// TreeR assumed to be be "prepared" in calling function
// by "MUON->GetTreeR(nev)" ????
TTree *treeR = fLoader->TreeR();
- fMUONData->SetTreeAddress("RC");
+ //AZ? fMUONData->SetTreeAddress("RC");
AddHitsForRecFromRawClusters(treeR);
// No sorting: it is done automatically in the previous function
}
TClonesArray *rawclusters;
AliDebug(1,"Enter AddHitsForRecFromRawClusters");
- nTRentries = Int_t(TR->GetEntries());
- if (nTRentries != 1) {
- AliError(Form("nTRentries = %d not equal to 1 ",nTRentries));
- exit(0);
+ if (fTrackMethod != 3) { //AZ
+ fMUONData->SetTreeAddress("RC"); //AZ
+ nTRentries = Int_t(TR->GetEntries());
+ if (nTRentries != 1) {
+ AliError(Form("nTRentries = %d not equal to 1 ",nTRentries));
+ exit(0);
+ }
+ fLoader->TreeR()->GetEvent(0); // only one entry
}
- fLoader->TreeR()->GetEvent(0); // only one entry
// Loop over tracking chambers
for (Int_t ch = 0; ch < AliMUONConstants::NTrackingCh(); ch++) {
(fNHitsForRecPerChamber[ch])++;
// more information into HitForRec
// resolution: info should be already in raw cluster and taken from it ????
- hitForRec->SetBendingReso2(fBendingResolution * fBendingResolution);
- hitForRec->SetNonBendingReso2(fNonBendingResolution * fNonBendingResolution);
+ //hitForRec->SetBendingReso2(fBendingResolution * fBendingResolution);
+ //hitForRec->SetNonBendingReso2(fNonBendingResolution * fNonBendingResolution);
+ hitForRec->SetBendingReso2(clus->GetErrY() * clus->GetErrY());
+ hitForRec->SetNonBendingReso2(clus->GetErrX() * clus->GetErrX());
// original raw cluster
hitForRec->SetChamberNumber(ch);
hitForRec->SetHitNumber(iclus);
hitForRec->Dump();}
} // end of cluster loop
} // end of chamber loop
- SortHitsForRecWithIncreasingChamber(); //AZ
+ SortHitsForRecWithIncreasingChamber();
return;
}
// To make the list of segments in all stations,
// from the list of hits to be reconstructed
AliDebug(1,"Enter MakeSegments");
- ResetSegments();
+ //AZ ResetSegments();
// Loop over stations
- Int_t nb = (fTrackMethod == 2) ? 3 : 0; //AZ
- //AZ for (Int_t st = 0; st < fgkMaxMuonTrackingStations; st++)
- for (Int_t st = nb; st < AliMUONConstants::NTrackingCh()/2; st++) //AZ
+ Int_t nb = (fTrackMethod != 1) ? 3 : 0; //AZ
+ for (Int_t st = nb; st < AliMUONConstants::NTrackingCh()/2; st++)
MakeSegmentsPerStation(st);
if (AliLog::GetGlobalDebugLevel() > 1) {
cout << "end of MakeSegments" << endl;
distNonBend = TMath::Abs(hit2Ptr->GetNonBendingCoor() - extNonBendCoor);
if (last2st) {
// bending slope
- bendingSlope = (hit1Ptr->GetBendingCoor() - hit2Ptr->GetBendingCoor()) /
- (hit1Ptr->GetZ() - hit2Ptr->GetZ());
- // absolute value of impact parameter
- impactParam =
- TMath::Abs(hit1Ptr->GetBendingCoor() - hit1Ptr->GetZ() * bendingSlope);
+ if ( hit1Ptr->GetZ() - hit2Ptr->GetZ() != 0.0 ) {
+ bendingSlope = (hit1Ptr->GetBendingCoor() - hit2Ptr->GetBendingCoor()) /
+ (hit1Ptr->GetZ() - hit2Ptr->GetZ());
+ // absolute value of impact parameter
+ impactParam =
+ TMath::Abs(hit1Ptr->GetBendingCoor() - hit1Ptr->GetZ() * bendingSlope);
+ }
+ else {
+ AliWarning("hit1Ptr->GetZ() = hit2Ptr->GetZ(): impactParam set to maxImpactParam");
+ impactParam = maxImpactParam;
+ }
}
// check for distances not too large,
// and impact parameter not too big if stations downstream of the dipole.
// from the list of segments and points in all stations
AliDebug(1,"Enter MakeTracks");
// The order may be important for the following Reset's
- ResetTracks();
- ResetTrackHits();
- if (fTrackMethod == 2) { //AZ - Kalman filter
+ //AZ ResetTracks();
+ //AZ ResetTrackHits();
+ if (fTrackMethod != 1) { //AZ - Kalman filter
MakeTrackCandidatesK();
if (fRecTracksPtr->GetEntriesFast() == 0) return;
// Follow tracks in stations(1..) 3, 2 and 1
TClonesArray *globalTrigger;
AliMUONLocalTrigger *locTrg;
AliMUONGlobalTrigger *gloTrg;
- AliMUONTriggerCircuit *circuit;
+
AliMUONTriggerTrack *recTriggerTrack = 0;
TTree* treeR = fLoader->TreeR();
// Loading MUON subsystem
AliMUON * pMUON = (AliMUON *) gAlice->GetDetector("MUON");
-
+ // Not really clean, but for the moment we must check whether the
+ // trigger uses new or old TriggerCircuit
+ Bool_t newTrigger=kFALSE;
+ if ( pMUON->DigitizerType().Contains("NewTrigger") ) newTrigger = kTRUE;
+
nTRentries = Int_t(treeR->GetEntries());
treeR->GetEvent(0); // only one entry
Float_t z11 = AliMUONConstants::DefaultChamberZ(10);
Float_t z21 = AliMUONConstants::DefaultChamberZ(12);
+ Float_t y11 = 0.;
+ Int_t stripX21 = 0;
+ Float_t y21 = 0.;
+ Float_t x11 = 0.;
+
for (Int_t i=0; i<nlocals; i++) { // loop on Local Trigger
locTrg = (AliMUONLocalTrigger*)localTrigger->UncheckedAt(i);
- circuit = &(pMUON->TriggerCircuit(locTrg->LoCircuit()));
- Float_t y11 = circuit->GetY11Pos(locTrg->LoStripX());
- Int_t stripX21 = locTrg->LoStripX()+locTrg->LoDev()+1;
- Float_t y21 = circuit->GetY21Pos(stripX21);
- Float_t x11 = circuit->GetX11Pos(locTrg->LoStripY());
+
+ if (!newTrigger) { // old trigger
+// printf("AliMUONTrackReconstructor::MakeTriggerTrack using OLD trigger \n");
+ AliMUONTriggerCircuit * circuit = &(pMUON->TriggerCircuit(locTrg->LoCircuit()));
+ y11 = circuit->GetY11Pos(locTrg->LoStripX());
+ stripX21 = locTrg->LoStripX()+locTrg->LoDev()+1;
+ y21 = circuit->GetY21Pos(stripX21);
+ x11 = circuit->GetX11Pos(locTrg->LoStripY());
+ } else { // new trigger
+// printf("AliMUONTrackReconstructor::MakeTriggerTrack using NEW trigger \n");
+ AliMUONTriggerCircuitNew * circuit =
+ &(pMUON->TriggerCircuitNew(locTrg->LoCircuit()-1)); // -1 !!!
+ y11 = circuit->GetY11Pos(locTrg->LoStripX());
+ stripX21 = locTrg->LoStripX()+locTrg->LoDev()+1;
+ y21 = circuit->GetY21Pos(stripX21);
+ x11 = circuit->GetX11Pos(locTrg->LoStripY());
+ }
+// printf(" MakeTriggerTrack %d %d %d %d %d %f %f %f \n",i,locTrg->LoCircuit(),locTrg->LoStripX(),locTrg->LoStripX()+locTrg->LoDev()+1,locTrg->LoStripY(),y11, y21, x11);
+
Float_t thetax = TMath::ATan2( x11 , z11 );
Float_t thetay = TMath::ATan2( (y21-y11) , (z21-z11) );
-
+
recTriggerTrack = new AliMUONTriggerTrack(x11,y11,thetax,thetay,gloTrigPat);
// since static statement does not work, set gloTrigPat for each track
-
+
fMUONData->AddRecTriggerTrack(*recTriggerTrack);
delete recTriggerTrack;
} // end of loop on Local Trigger
Int_t ch = -1, chInStation, chBestHit = -1, iHit, iSegment, station, trackIndex;
Double_t bestChi2, chi2, dZ1, dZ2, dZ3, maxSigma2Distance, mcsFactor;
Double_t bendingMomentum, chi2Norm = 0.;
- AliMUON *pMUON = (AliMUON*) gAlice->GetModule("MUON"); // necessary ????
+
// local maxSigma2Distance, for easy increase in testing
maxSigma2Distance = fMaxSigma2Distance;
AliDebug(2,"Enter FollowTracks");
mcsFactor = 0.0136 * trackParam1->GetInverseBendingMomentum();
mcsFactor = fChamberThicknessInX0 * mcsFactor * mcsFactor;
// Z difference from previous station
- dZ1 = (&(pMUON->Chamber(2 * station)))->Z() -
- (&(pMUON->Chamber(2 * station + 2)))->Z();
+ dZ1 = AliMUONConstants::DefaultChamberZ(2 * station) -
+ AliMUONConstants::DefaultChamberZ(2 * station + 2);
// Z difference between the two previous stations
- dZ2 = (&(pMUON->Chamber(2 * station + 2)))->Z() -
- (&(pMUON->Chamber(2 * station + 4)))->Z();
+ dZ2 = AliMUONConstants::DefaultChamberZ(2 * station + 2) -
+ AliMUONConstants::DefaultChamberZ(2 * station + 4);
// Z difference between the two chambers in the previous station
- dZ3 = (&(pMUON->Chamber(2 * station)))->Z() -
- (&(pMUON->Chamber(2 * station + 1)))->Z();
+ dZ3 = AliMUONConstants::DefaultChamberZ(2 * station) -
+ AliMUONConstants::DefaultChamberZ(2 * station + 1);
extrapSegment->SetBendingCoorReso2(fBendingResolution * fBendingResolution);
extrapSegment->
SetNonBendingCoorReso2(fNonBendingResolution * fNonBendingResolution);
AliMUONTrackK *trackK;
AliDebug(1,"Enter MakeTrackCandidatesK");
- // Reset the TClonesArray of reconstructed tracks
- if (fRecTracksPtr) fRecTracksPtr->Delete();
- // Delete in order that the Track destructors are called,
- // hence the space for the TClonesArray of pointers to TrackHit's is freed
- fNRecTracks = 0;
- AliMUONTrackK a(this, fHitsForRecPtr); // bad idea ???
+ AliMUONTrackK a(this, fHitsForRecPtr);
// Loop over stations(1...) 5 and 4
for (istat=4; istat>=3; istat--) {
// Loop over segments in the station
for (iseg=0; iseg<fNSegments[istat]; iseg++) {
// Transform segments to tracks and evaluate covariance matrix
segment = (AliMUONSegment*) ((*fSegmentsPtr[istat])[iseg]);
- trackK = new ((*fRecTracksPtr)[fNRecTracks]) AliMUONTrackK(segment);
- fNRecTracks++;
+ trackK = new ((*fRecTracksPtr)[fNRecTracks++]) AliMUONTrackK(segment);
} // for (iseg=0;...)
} // for (istat=4;...)
return;
TClonesArray *rawclusters;
clus = 0; rawclusters = 0;
- //AZ(z->-z) zDipole1 = GetSimpleBPosition() - GetSimpleBLength()/2;
- //AZ(z->-z) zDipole2 = zDipole1 + GetSimpleBLength();
zDipole1 = GetSimpleBPosition() + GetSimpleBLength()/2;
zDipole2 = zDipole1 - GetSimpleBLength();
rawclusters = fMUONData->RawClusters(hit->GetChamberNumber());
clus = (AliMUONRawCluster*) rawclusters->UncheckedAt(hit->
GetHitNumber());
- printf("%3d", clus->GetTrack(1)-1);
- if (clus->GetTrack(2) != 0) printf("%3d \n", clus->GetTrack(2)-1);
+ printf(" %d", clus->GetTrack(1));
+ if (clus->GetTrack(2) != -1) printf(" %d \n", clus->GetTrack(2));
else printf("\n");
} // if (fRecTrackRefHits)
}
ok = kTRUE;
if (trackK->GetRecover() == 0) hit = (AliMUONHitForRec*)
- trackK->GetHitOnTrack()->Last(); // last hit
- //else hit = (AliMUONHitForRec*) (*trackK->GetHitOnTrack())[1]; // 2'nd hit
+ trackK->GetTrackHits()->Last(); // last hit
else hit = trackK->GetHitLastOk(); // hit where track stopped
if (hit) ichamBeg = hit->GetChamberNumber();
ichamEnd = 0;
// Check propagation direction
- if (hit == NULL) ichamBeg = ichamEnd;
- //AZ(z->-z) else if (trackK->GetTrackDir() > 0) {
+ if (!hit) { ichamBeg = ichamEnd; AliFatal(" ??? "); }
else if (trackK->GetTrackDir() < 0) {
ichamEnd = 9; // forward propagation
ok = trackK->KalmanFilter(ichamBeg,ichamEnd,kFALSE,zDipole1,zDipole2);
ichamEnd = 6; // backward propagation
// Change weight matrix and zero fChi2 for backpropagation
trackK->StartBack();
- //AZ trackK->SetTrackDir(-1);
trackK->SetTrackDir(1);
ok = trackK->KalmanFilter(ichamBeg,ichamEnd,kTRUE,zDipole1,zDipole2);
ichamBeg = ichamEnd;
}
if (ok) {
- //AZ trackK->SetTrackDir(-1);
trackK->SetTrackDir(1);
trackK->SetBPFlag(kFALSE);
ok = trackK->KalmanFilter(ichamBeg,ichamEnd,kFALSE,zDipole1,zDipole2);
chamBits = 0;
Double_t chi2max = 0;
for (Int_t i=0; i<trackK->GetNTrackHits(); i++) {
- hit = (AliMUONHitForRec*) (*trackK->GetHitOnTrack())[i];
+ hit = (AliMUONHitForRec*) (*trackK->GetTrackHits())[i];
chamBits |= BIT(hit->GetChamberNumber());
if (trackK->GetChi2PerPoint(i) > chi2max) chi2max = trackK->GetChi2PerPoint(i);
}
AliMUONHitForRec *hit1, *hit2, *hit;
track1 = (AliMUONTrackK*) ((*fRecTracksPtr)[icand]);
- hit1 = (AliMUONHitForRec*) (*track1->GetHitOnTrack())[0]; // 1'st hit
- hit2 = (AliMUONHitForRec*) (*track1->GetHitOnTrack())[1]; // 2'nd hit
+ hit1 = (AliMUONHitForRec*) (*track1->GetTrackHits())[0]; // 1'st hit
+ hit2 = (AliMUONHitForRec*) (*track1->GetTrackHits())[1]; // 2'nd hit
for (Int_t i=0; i<icand; i++) {
track2 = (AliMUONTrackK*) ((*fRecTracksPtr)[i]);
} else {
Int_t nSame = 0;
for (Int_t j=0; j<track2->GetNTrackHits(); j++) {
- hit = (AliMUONHitForRec*) (*track2->GetHitOnTrack())[j];
+ hit = (AliMUONHitForRec*) (*track2->GetTrackHits())[j];
if (hit == hit1 || hit == hit2) {
nSame++;
if (nSame == 2) return kFALSE;
{
// Set track method and recreate track container if necessary
- fTrackMethod = iTrackMethod == 2 ? 2 : 1;
- if (fTrackMethod == 2) {
+ fTrackMethod = TMath::Min (iTrackMethod, 3);
+ fTrackMethod = TMath::Max (fTrackMethod, 1);
+ if (fTrackMethod != 1) {
if (fRecTracksPtr) delete fRecTracksPtr;
fRecTracksPtr = new TClonesArray("AliMUONTrackK", 10);
- cout << " *** Tracking with the Kalman filter *** " << endl;
- } else cout << " *** Traditional tracking *** " << endl;
+ if (fTrackMethod == 2) cout << " *** Tracking with the Kalman filter *** " << endl;
+ else cout << " *** Combined cluster / track finder ***" << endl;
+ } else cout << " *** Original tracking *** " << endl;
}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff