X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=MUON%2FAliMUONRecoCheck.cxx;h=144d0769280a9f795fa56cbca853bdc216bc0ac6;hb=542e9a388e60e7b3a7ee53bc3ad4118c0463c8d1;hp=b1b800d6a9355b5bf51bbe5c7f7e87a988f957e6;hpb=f29ba3e1555f9cc21a4e8926e589d23342306213;p=u%2Fmrichter%2FAliRoot.git diff --git a/MUON/AliMUONRecoCheck.cxx b/MUON/AliMUONRecoCheck.cxx index b1b800d6a93..144d0769280 100644 --- a/MUON/AliMUONRecoCheck.cxx +++ b/MUON/AliMUONRecoCheck.cxx @@ -1,465 +1,471 @@ /************************************************************************** - * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * - * * - * Author: The ALICE Off-line Project. * - * Contributors are mentioned in the code where appropriate. * - * * - * Permission to use, copy, modify and distribute this software and its * - * documentation strictly for non-commercial purposes is hereby granted * - * without fee, provided that the above copyright notice appears in all * - * copies and that both the copyright notice and this permission notice * - * appear in the supporting documentation. The authors make no claims * - * about the suitability of this software for any purpose. It is * - * provided "as is" without express or implied warranty. * - **************************************************************************/ +* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * +* * +* Author: The ALICE Off-line Project. * +* Contributors are mentioned in the code where appropriate. * +* * +* Permission to use, copy, modify and distribute this software and its * +* documentation strictly for non-commercial purposes is hereby granted * +* without fee, provided that the above copyright notice appears in all * +* copies and that both the copyright notice and this permission notice * +* appear in the supporting documentation. The authors make no claims * +* about the suitability of this software for any purpose. It is * +* provided "as is" without express or implied warranty. * +**************************************************************************/ + +/* $Id$ */ + +//----------------------------------------------------------------------------- +/// \class AliMUONRecoCheck +/// Utility class to check reconstruction +/// Reconstructed tracks are compared to reference tracks. +/// The reference tracks are built from AliTrackReference for the +/// hit in chamber (0..9) and from kinematics for the vertex parameters. +//----------------------------------------------------------------------------- + +// 13 Nov 2007: +// Added a method to create a list of reconstructed AliMUONTrack objects from +// ESD data. This is necessary since the track objects that are actually created +// during offline reconstruction are no longer stored to disk. +// - Artur Szostak -////////////////////////////////////////////////////////////////////////////////////// -// -// Utility class to check the muon reconstruction. Reconstructed tracks are compared -// to reference tracks. The reference tracks are built from AliTrackReference for the -// hit in chamber (0..9) and from kinematics for the vertex parameters. -// -////////////////////////////////////////////////////////////////////////////////////// - - -#include - -#include "AliRun.h" // for gAlice -#include "AliLog.h" -#include "AliLoader.h" -#include "AliRunLoader.h" -#include "AliTrackReference.h" -#include "AliHeader.h" -#include "AliMC.h" -#include "AliStack.h" -#include "AliMUON.h" #include "AliMUONRecoCheck.h" +#include "AliMUONRawClusterV2.h" #include "AliMUONTrack.h" -#include "AliMUONData.h" -#include "AliMUONChamber.h" #include "AliMUONConstants.h" +#include "AliMUONTrackStoreV1.h" + +#include "AliMCEventHandler.h" +#include "AliMCEvent.h" +#include "AliStack.h" +#include "AliTrackReference.h" +#include "AliLog.h" +#include "AliESDEvent.h" +#include "AliESDMuonTrack.h" +#include +#include +#include +#include +#include + +/// \cond CLASSIMP ClassImp(AliMUONRecoCheck) +/// \endcond //_____________________________________________________________________________ -AliMUONRecoCheck::AliMUONRecoCheck(Char_t *chLoader) +AliMUONRecoCheck::AliMUONRecoCheck(Char_t *esdFileName, Char_t *pathSim) +: TObject(), +fMCEventHandler(new AliMCEventHandler()), +fESDEvent(new AliESDEvent()), +fESDTree (0x0), +fESDFile (0x0), +fCurrentEvent(0), +fTrackRefStore(0x0), +fRecoTrackRefStore(0x0), +fRecoTrackStore(0x0) { -// Constructor - fMuonTrackRef = new TClonesArray("AliMUONTrack", 10); - - // open the run loader - fRunLoader = AliRunLoader::Open(chLoader); - if (!fRunLoader) { - AliError(Form("no run loader found in file %s","galice.root" )); + /// Normal ctor + + // TrackRefs and Particules + fMCEventHandler->SetInputPath(pathSim); + fMCEventHandler->InitIO(""); + + // ESD MUON Tracks + fESDFile = TFile::Open(esdFileName); // open the file + if (!fESDFile || !fESDFile->IsOpen()) { + AliError(Form("opening ESD file %s failed", esdFileName)); + fESDFile = 0x0; return; } - // initialize loader's - AliLoader *loader = fRunLoader->GetLoader("MUONLoader"); - - // initialize container - fMUONData = new AliMUONData(loader,"MUON","MUON"); - - // Loading AliRun master - if (fRunLoader->GetAliRun() == 0x0) fRunLoader->LoadgAlice(); - - fRunLoader->LoadKinematics("READ"); - fRunLoader->LoadTrackRefs("READ"); - loader->LoadTracks("READ"); + fESDTree = (TTree*) fESDFile->Get("esdTree"); // get the tree + if (!fESDTree) { + AliError("no ESD tree found"); + fESDFile->Close(); + fESDFile = 0x0; + return; + } + fESDEvent->ReadFromTree(fESDTree); // link fESDEvent to the tree +} - fReconstructibleTracks = 0; - fRecoTracks = 0; +//_____________________________________________________________________________ +AliMUONRecoCheck::~AliMUONRecoCheck() +{ + /// Destructor + delete fMCEventHandler; + delete fESDEvent; + if (fESDFile) fESDFile->Close(); + ResetStores(); } -//____________________________________________________________________ -AliMUONRecoCheck::AliMUONRecoCheck(const AliMUONRecoCheck& rhs) - : TObject(rhs) +//_____________________________________________________________________________ +void AliMUONRecoCheck::ResetStores() { -// Protected copy constructor + /// Deletes all the store objects that have been created and resets the pointers to 0x0 + delete fTrackRefStore; fTrackRefStore = 0x0; + delete fRecoTrackRefStore; fRecoTrackRefStore = 0x0; + delete fRecoTrackStore; fRecoTrackStore = 0x0; +} - AliFatal("Not implemented."); +//_____________________________________________________________________________ +Int_t AliMUONRecoCheck::NumberOfEvents() const +{ + /// Return the number of events + if (fESDTree) return fESDTree->GetEntries(); + return 0; } //_____________________________________________________________________________ -AliMUONRecoCheck::~AliMUONRecoCheck() +AliMUONVTrackStore* AliMUONRecoCheck::ReconstructedTracks(Int_t event) { - fRunLoader->UnloadKinematics(); - fRunLoader->UnloadTrackRefs(); - fRunLoader->UnloadTracks(); - fMuonTrackRef->Delete(); - delete fMuonTrackRef; - delete fMUONData; + /// Return a track store containing the reconstructed tracks (converted into + /// MUONTrack objects) for a given event + if (event != fCurrentEvent) { + ResetStores(); + fCurrentEvent = event; + } + + if (fRecoTrackStore != 0x0) return fRecoTrackStore; + else { + if (!fESDTree) return 0x0; + if (fESDTree->GetEvent(event) <= 0) { + AliError(Form("fails to read ESD object for event %d", event)); + return 0x0; + } + MakeReconstructedTracks(); + return fRecoTrackStore; + } } -//________________________________________________________________________ -AliMUONRecoCheck& AliMUONRecoCheck::operator = (const AliMUONRecoCheck& rhs) +//_____________________________________________________________________________ +AliMUONVTrackStore* AliMUONRecoCheck::TrackRefs(Int_t event) { -// Protected assignement operator + /// Return a track store containing the track references (converted into + /// MUONTrack objects) for a given event + if (event != fCurrentEvent) { + ResetStores(); + fCurrentEvent = event; + } + + if (fTrackRefStore != 0x0) return fTrackRefStore; + else { + if (!fMCEventHandler->GetEvent(event)) { + AliError(Form("fails to read MC objects for event %d", event)); + return 0x0; + } + MakeTrackRefs(); + return fTrackRefStore; + } +} - if (this == &rhs) return *this; +//_____________________________________________________________________________ +AliMUONVTrackStore* AliMUONRecoCheck::ReconstructibleTracks(Int_t event) +{ + /// Return a track store containing the reconstructible tracks for a given event + if (event != fCurrentEvent) { + ResetStores(); + fCurrentEvent = event; + } + + if (fRecoTrackRefStore != 0x0) return fRecoTrackRefStore; + else { + if (TrackRefs(event) == 0x0) return 0x0; + MakeReconstructibleTracks(); + return fRecoTrackRefStore; + } +} - AliFatal("Not implemented."); +//_____________________________________________________________________________ +void AliMUONRecoCheck::MakeReconstructedTracks() +{ + /// Make reconstructed tracks + fRecoTrackStore = new AliMUONTrackStoreV1(); + + Int_t nTracks = Int_t(fESDEvent->GetNumberOfMuonTracks()); + + // loop over all reconstructed tracks + for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) { + AliESDMuonTrack* esdTrack = fESDEvent->GetMuonTrack(iTrack); + + // copy ESD MUON track into standard MUON track + AliMUONTrack track(*esdTrack); - return *this; + // add track to the store + fRecoTrackStore->Add(track); + } + } //_____________________________________________________________________________ -void AliMUONRecoCheck::MakeTrackRef() +void AliMUONRecoCheck::MakeTrackRefs() { - // Make reconstructible tracks - - AliTrackReference *trackReference; - AliMUONTrack *muonTrack; - AliMUONTrackParam *trackParam; - AliMUONHitForRec *hitForRec; - Float_t x, y, z, pX, pY, pZ, pYZ; - Int_t track, trackSave; - TParticle *particle; - Float_t bendingSlope = 0; - Float_t nonBendingSlope = 0; - Float_t inverseBendingMomentum = 0; + /// Make reconstructible tracks + AliMUONVTrackStore *tmpTrackRefStore = new AliMUONTrackStoreV1(); - TTree* treeTR = fRunLoader->TreeTR(); - if (treeTR == NULL) return; - - TBranch* branch = treeTR->GetBranch("MUON"); - if (branch == NULL) return; - - TClonesArray* trackRefs = 0; - branch->SetAddress(&trackRefs); - branch->SetAutoDelete(kTRUE); - - Int_t nTrackRef = (Int_t)branch->GetEntries(); - - track = trackSave = -999; - Bool_t isNewTrack; - Int_t iHitMin, iChamber; - - trackParam = new AliMUONTrackParam(); - hitForRec = new AliMUONHitForRec(); - muonTrack = new AliMUONTrack(); - - for (Int_t iTrackRef = 0; iTrackRef < nTrackRef; iTrackRef++) { - branch->GetEntry(iTrackRef); + Double_t x, y, z, pX, pY, pZ, bendingSlope, nonBendingSlope, inverseBendingMomentum; + TParticle* particle; + TClonesArray* trackRefs; + Int_t nTrackRef = fMCEventHandler->MCEvent()->GetNumberOfTracks(); + + // loop over simulated tracks + for (Int_t iTrackRef = 0; iTrackRef < nTrackRef; ++iTrackRef) { + Int_t nHits = fMCEventHandler->GetParticleAndTR(iTrackRef, particle, trackRefs); - iHitMin = 0; - isNewTrack = kTRUE; + // skip empty trackRefs + if (nHits < 1) continue; - if (!trackRefs->GetEntries()) continue; - - while (isNewTrack) { + // get the particle charge for further calculation + TParticlePDG* ppdg = particle->GetPDG(); + Int_t charge = (Int_t)(ppdg->Charge()/3.0); + + AliMUONTrack track; + + // loop over simulated track hits + for (Int_t iHit = 0; iHit < nHits; ++iHit) { + AliTrackReference* trackReference = static_cast(trackRefs->UncheckedAt(iHit)); - for (Int_t iHit = iHitMin; iHit < trackRefs->GetEntries(); iHit++) { + // skip trackRefs not in MUON + if (trackReference->DetectorId() != AliTrackReference::kMUON) continue; + + // Get track parameters of current hit + x = trackReference->X(); + y = trackReference->Y(); + z = trackReference->Z(); + pX = trackReference->Px(); + pY = trackReference->Py(); + pZ = trackReference->Pz(); - trackReference = (AliTrackReference*)trackRefs->At(iHit); - x = trackReference->X(); - y = trackReference->Y(); - z = trackReference->Z(); - pX = trackReference->Px(); - pY = trackReference->Py(); - pZ = trackReference->Pz(); - - track = trackReference->GetTrack(); - - if (track != trackSave && iHit != 0) { - iHitMin = iHit; - trackSave = track; - break; - } - - // track parameters at hit - trackParam->SetBendingCoor(y); - trackParam->SetNonBendingCoor(x); - trackParam->SetZ(z); - - if (TMath::Abs(pZ) > 0) { - bendingSlope = pY/pZ; - nonBendingSlope = pX/pZ; - } - pYZ = TMath::Sqrt(pY*pY+pZ*pZ); - if (pYZ >0) inverseBendingMomentum = 1/pYZ; - - trackParam->SetBendingSlope(bendingSlope); - trackParam->SetNonBendingSlope(nonBendingSlope); - trackParam->SetInverseBendingMomentum(inverseBendingMomentum); - - hitForRec->SetBendingCoor(y); - hitForRec->SetNonBendingCoor(x); - hitForRec->SetZ(z); - hitForRec->SetBendingReso2(0.0); - hitForRec->SetNonBendingReso2(0.0); - iChamber = AliMUONConstants::ChamberNumber(z); - hitForRec->SetChamberNumber(iChamber); - - muonTrack->AddTrackParamAtHit(trackParam); - muonTrack->AddHitForRecAtHit(hitForRec); - muonTrack->SetTrackID(track); - - trackSave = track; - if (iHit == trackRefs->GetEntries()-1) isNewTrack = kFALSE; - } + // check chamberId of current trackReference + Int_t chamberId = AliMUONConstants::ChamberNumber(z); + if (chamberId < 0 || chamberId >= AliMUONConstants::NTrackingCh()) continue; - // track parameters at vertex - particle = fRunLoader->GetHeader()->Stack()->Particle(muonTrack->GetTrackID()); - if (particle) { - x = particle->Vx(); - y = particle->Vy(); - z = particle->Vz(); - pX = particle->Px(); - pY = particle->Py(); - pZ = particle->Pz(); - - trackParam->SetBendingCoor(y); - trackParam->SetNonBendingCoor(x); - trackParam->SetZ(z); - if (TMath::Abs(pZ) > 0) { - bendingSlope = pY/pZ; - nonBendingSlope = pX/pZ; - } - pYZ = TMath::Sqrt(pY*pY+pZ*pZ); - if (pYZ >0) inverseBendingMomentum = 1/pYZ; - trackParam->SetBendingSlope(bendingSlope); - trackParam->SetNonBendingSlope(nonBendingSlope); - trackParam->SetInverseBendingMomentum(inverseBendingMomentum); + // set hit parameters + AliMUONRawClusterV2 hit(chamberId, 0, 0); + hit.SetXYZ(x,y,z); + hit.SetErrXY(0.,0.); - muonTrack->SetTrackParamAtVertex(trackParam); + // compute track parameters at hit + Double_t bendingSlope = 0; + Double_t nonBendingSlope = 0; + Double_t inverseBendingMomentum = 0; + if (TMath::Abs(pZ) > 0) { + bendingSlope = pY/pZ; + nonBendingSlope = pX/pZ; } + Double_t pYZ = TMath::Sqrt(pY*pY+pZ*pZ); + if (pYZ >0) inverseBendingMomentum = 1/pYZ; + inverseBendingMomentum *= charge; - AddMuonTrackReference(muonTrack); - muonTrack->ResetTrackParamAtHit(); - muonTrack->ResetHitForRecAtHit(); + // set track parameters at hit + AliMUONTrackParam trackParam; + trackParam.SetNonBendingCoor(x); + trackParam.SetBendingCoor(y); + trackParam.SetZ(z); + trackParam.SetBendingSlope(bendingSlope); + trackParam.SetNonBendingSlope(nonBendingSlope); + trackParam.SetInverseBendingMomentum(inverseBendingMomentum); - } // end while isNewTrack + // add track parameters at current hit to the track + track.AddTrackParamAtCluster(trackParam,hit,kTRUE); + } + + // if none of the track hits was in MUON, goto the next track + if (track.GetNClusters() < 1) continue; + + // get track parameters at particle's vertex + x = particle->Vx(); + y = particle->Vy(); + z = particle->Vz(); + pX = particle->Px(); + pY = particle->Py(); + pZ = particle->Pz(); + // compute rest of track parameters at particle's vertex + bendingSlope = 0; + nonBendingSlope = 0; + inverseBendingMomentum = 0; + if (TMath::Abs(pZ) > 0) { + bendingSlope = pY/pZ; + nonBendingSlope = pX/pZ; + } + Double_t pYZ = TMath::Sqrt(pY*pY+pZ*pZ); + if (pYZ >0) inverseBendingMomentum = 1/pYZ; + inverseBendingMomentum *= charge; + + // set track parameters at particle's vertex + AliMUONTrackParam trackParamAtVertex; + trackParamAtVertex.SetNonBendingCoor(x); + trackParamAtVertex.SetBendingCoor(y); + trackParamAtVertex.SetZ(z); + trackParamAtVertex.SetBendingSlope(bendingSlope); + trackParamAtVertex.SetNonBendingSlope(nonBendingSlope); + trackParamAtVertex.SetInverseBendingMomentum(inverseBendingMomentum); + + // add track parameters at vertex + track.SetTrackParamAtVertex(&trackParamAtVertex); + + // sort trackParamAtCluster and store the track + track.GetTrackParamAtCluster()->Sort(); + track.SetTrackID(iTrackRef); + tmpTrackRefStore->Add(track); } - CleanMuonTrackRef(); + CleanMuonTrackRef(tmpTrackRefStore); - ReconstructibleTracks(); - - delete muonTrack; - delete trackParam; - delete hitForRec; - + delete tmpTrackRefStore; } -//____________________________________________________________________________ -TClonesArray* AliMUONRecoCheck::GetTrackReco() -{ - // Return TClonesArray of reconstructed tracks - - GetMUONData()->ResetRecTracks(); - GetMUONData()->SetTreeAddress("RT"); - fTrackReco = GetMUONData()->RecTracks(); - GetMUONData()->GetRecTracks(); - fRecoTracks = fTrackReco->GetEntriesFast(); - return fTrackReco; -} //_____________________________________________________________________________ -void AliMUONRecoCheck::PrintEvent() const +void AliMUONRecoCheck::CleanMuonTrackRef(const AliMUONVTrackStore *tmpTrackRefStore) { - // debug facility - - AliMUONTrack *track; - AliMUONHitForRec *hitForRec; - TClonesArray * hitForRecAtHit = 0; - Float_t xRec,yRec,zRec; - - Int_t nTrackRef = fMuonTrackRef->GetEntriesFast(); + /// Re-calculate hits parameters because two AliTrackReferences are recorded for + /// each chamber (one when particle is entering + one when particle is leaving + /// the sensitive volume) + fTrackRefStore = new AliMUONTrackStoreV1(); - printf(" ******************************* \n"); - printf(" nb of tracks %d \n",nTrackRef); - - for (Int_t index = 0; index < nTrackRef; index++) { - track = (AliMUONTrack*)fMuonTrackRef->At(index); - hitForRecAtHit = track->GetHitForRecAtHit(); - Int_t nTrackHits = hitForRecAtHit->GetEntriesFast(); - printf(" track number %d \n",index); - for (Int_t iHit = 0; iHit < nTrackHits; iHit++){ - hitForRec = (AliMUONHitForRec*) hitForRecAtHit->At(iHit); - xRec = hitForRec->GetNonBendingCoor(); - yRec = hitForRec->GetBendingCoor(); - zRec = hitForRec->GetZ(); - printf(" x,y,z: %f , %f , %f \n",xRec,yRec,zRec); - } - } -} - -//_____________________________________________________________________________ -void AliMUONRecoCheck::ResetTracks() const -{ - if (fMuonTrackRef) fMuonTrackRef->Clear(); -} -//_____________________________________________________________________________ -void AliMUONRecoCheck::CleanMuonTrackRef() -{ - // Re-calculate hits parameters because two AliTrackReferences are recorded for - // each chamber (one when particle is entering + one when particle is leaving - // the sensitive volume) - - Float_t maxGasGap = 1.; // cm - AliMUONTrack *track, *trackNew; - AliMUONHitForRec *hitForRec, *hitForRec1, *hitForRec2; - AliMUONTrackParam *trackParam, *trackParam1, *trackParam2, *trackParamAtVertex; - TClonesArray * hitForRecAtHit = 0; - TClonesArray * trackParamAtHit = 0; - Float_t xRec,yRec,zRec; - Float_t xRec1,yRec1,zRec1; - Float_t xRec2,yRec2,zRec2; - Float_t bendingSlope,nonBendingSlope,bendingMomentum; - Float_t bendingSlope1,nonBendingSlope1,bendingMomentum1; - Float_t bendingSlope2,nonBendingSlope2,bendingMomentum2; - TClonesArray *newMuonTrackRef = new TClonesArray("AliMUONTrack", 10); - Int_t iHit1; - Int_t iChamber = 0; - Int_t nRec = 0; - Int_t nTrackHits = 0; - - hitForRec = new AliMUONHitForRec(); - trackParam = new AliMUONTrackParam(); - trackNew = new AliMUONTrack(); - - Int_t nTrackRef = fMuonTrackRef->GetEntriesFast(); + Double_t maxGasGap = 1.; // cm + Double_t x, y, z, pX, pY, pZ, x1, y1, z1, pX1, pY1, pZ1, z2; + Double_t bendingSlope,nonBendingSlope,inverseBendingMomentum; - for (Int_t index = 0; index < nTrackRef; index++) { - track = (AliMUONTrack*)fMuonTrackRef->At(index); - hitForRecAtHit = track->GetHitForRecAtHit(); - trackParamAtHit = track->GetTrackParamAtHit(); - trackParamAtVertex = track->GetTrackParamAtVertex(); - nTrackHits = hitForRecAtHit->GetEntriesFast(); - iHit1 = 0; + // create iterator + TIter next(tmpTrackRefStore->CreateIterator()); + + // loop over tmpTrackRef + AliMUONTrack* track; + while ( ( track = static_cast(next()) ) ) { + + AliMUONTrack newTrack; + + // loop over tmpTrackRef's hits + Int_t iHit1 = 0; + Int_t nTrackHits = track->GetNClusters(); while (iHit1 < nTrackHits) { - hitForRec1 = (AliMUONHitForRec*) hitForRecAtHit->At(iHit1); - trackParam1 = (AliMUONTrackParam*) trackParamAtHit->At(iHit1); - xRec1 = hitForRec1->GetNonBendingCoor(); - yRec1 = hitForRec1->GetBendingCoor(); - zRec1 = hitForRec1->GetZ(); - xRec = xRec1; - yRec = yRec1; - zRec = zRec1; - bendingSlope1 = trackParam1->GetBendingSlope(); - nonBendingSlope1 = trackParam1->GetNonBendingSlope(); - bendingMomentum1 = 0; - if (TMath::Abs(trackParam1->GetInverseBendingMomentum()) > 0) - bendingMomentum1 = 1./trackParam1->GetInverseBendingMomentum(); - bendingSlope = bendingSlope1; - nonBendingSlope = nonBendingSlope1; - bendingMomentum = bendingMomentum1; - nRec = 1; + AliMUONTrackParam *trackParam1 = (AliMUONTrackParam*) track->GetTrackParamAtCluster()->UncheckedAt(iHit1); + + // get track parameters at hit1 + x1 = trackParam1->GetNonBendingCoor(); + y1 = trackParam1->GetBendingCoor(); + z1 = trackParam1->GetZ(); + pX1 = trackParam1->Px(); + pY1 = trackParam1->Py(); + pZ1 = trackParam1->Pz(); + + // prepare new track parameters + x = x1; + y = y1; + z = z1; + pX = pX1; + pY = pY1; + pZ = pZ1; + + // loop over next tmpTrackRef's hits + Int_t nCombinedHits = 1; for (Int_t iHit2 = iHit1+1; iHit2 < nTrackHits; iHit2++) { - hitForRec2 = (AliMUONHitForRec*) hitForRecAtHit->At(iHit2); - trackParam2 = (AliMUONTrackParam*) trackParamAtHit->At(iHit2); - xRec2 = hitForRec2->GetNonBendingCoor(); - yRec2 = hitForRec2->GetBendingCoor(); - zRec2 = hitForRec2->GetZ(); - bendingSlope2 = trackParam2->GetBendingSlope(); - nonBendingSlope2 = trackParam2->GetNonBendingSlope(); - bendingMomentum2 = 0; - if (TMath::Abs(trackParam2->GetInverseBendingMomentum()) > 0) - bendingMomentum2 = 1./trackParam2->GetInverseBendingMomentum(); - - if ( TMath::Abs(zRec2-zRec1) < maxGasGap ) { - nRec++; - xRec += xRec2; - yRec += yRec2; - zRec += zRec2; - bendingSlope += bendingSlope2; - nonBendingSlope += nonBendingSlope2; - bendingMomentum += bendingMomentum2; - iHit1 = iHit2; - } + AliMUONTrackParam *trackParam2 = (AliMUONTrackParam*) track->GetTrackParamAtCluster()->UncheckedAt(iHit2); + + // get z position of hit2 + z2 = trackParam2->GetZ(); - } // end iHit2 - xRec /= (Float_t)nRec; - yRec /= (Float_t)nRec; - zRec /= (Float_t)nRec; - bendingSlope /= (Float_t)nRec; - nonBendingSlope /= (Float_t)nRec; - bendingMomentum /= (Float_t)nRec; + // complete new track parameters if hit2 is on the same detection element + if ( TMath::Abs(z2-z1) < maxGasGap ) { + x += trackParam2->GetNonBendingCoor(); + y += trackParam2->GetBendingCoor(); + z += z2; + pX += trackParam2->Px(); + pY += trackParam2->Py(); + pZ += trackParam2->Pz(); + nCombinedHits++; + iHit1 = iHit2; + } + + } - hitForRec->SetNonBendingCoor(xRec); - hitForRec->SetBendingCoor(yRec); - hitForRec->SetZ(zRec); - iChamber = AliMUONConstants::ChamberNumber(zRec); - hitForRec->SetChamberNumber(iChamber); - hitForRec->SetBendingReso2(0.0); - hitForRec->SetNonBendingReso2(0.0); - trackParam->SetNonBendingCoor(xRec); - trackParam->SetBendingCoor(yRec); - trackParam->SetZ(zRec); - trackParam->SetNonBendingSlope(nonBendingSlope); - trackParam->SetBendingSlope(bendingSlope); - if (TMath::Abs(bendingMomentum) > 0) - trackParam->SetInverseBendingMomentum(1./bendingMomentum); - - trackNew->AddHitForRecAtHit(hitForRec); - trackNew->AddTrackParamAtHit(trackParam); + // finalize new track parameters + x /= (Double_t)nCombinedHits; + y /= (Double_t)nCombinedHits; + z /= (Double_t)nCombinedHits; + pX /= (Double_t)nCombinedHits; + pY /= (Double_t)nCombinedHits; + pZ /= (Double_t)nCombinedHits; + bendingSlope = 0; + nonBendingSlope = 0; + inverseBendingMomentum = 0; + if (TMath::Abs(pZ) > 0) { + bendingSlope = pY/pZ; + nonBendingSlope = pX/pZ; + } + Double_t pYZ = TMath::Sqrt(pY*pY+pZ*pZ); + if (pYZ >0) inverseBendingMomentum = 1/pYZ; + inverseBendingMomentum *= trackParam1->GetCharge(); + + // set hit parameters + AliMUONRawClusterV2 hit(trackParam1->GetClusterPtr()->GetChamberId(), 0, 0); + hit.SetXYZ(x,y,z); + hit.SetErrXY(0.,0.); + + // set new track parameters at new hit + AliMUONTrackParam trackParam; + trackParam.SetNonBendingCoor(x); + trackParam.SetBendingCoor(y); + trackParam.SetZ(z); + trackParam.SetBendingSlope(bendingSlope); + trackParam.SetNonBendingSlope(nonBendingSlope); + trackParam.SetInverseBendingMomentum(inverseBendingMomentum); + + // add track parameters at current hit to the track + newTrack.AddTrackParamAtCluster(trackParam,hit,kTRUE); iHit1++; - } // end iHit1 - - trackNew->SetTrackID(track->GetTrackID()); - trackNew->SetTrackParamAtVertex(trackParamAtVertex); - {new ((*newMuonTrackRef)[newMuonTrackRef->GetEntriesFast()]) AliMUONTrack(*trackNew);} - trackNew->ResetHitForRecAtHit(); - trackNew->ResetTrackParamAtHit(); + } + + newTrack.GetTrackParamAtCluster()->Sort(); + newTrack.SetTrackID(track->GetTrackID()); + newTrack.SetTrackParamAtVertex(track->GetTrackParamAtVertex()); + fTrackRefStore->Add(newTrack); - } // end trackRef - - fMuonTrackRef->Clear(); - nTrackRef = newMuonTrackRef->GetEntriesFast(); - for (Int_t index = 0; index < nTrackRef; index++) { - track = (AliMUONTrack*)newMuonTrackRef->At(index); - AddMuonTrackReference(track); } - - delete trackNew; - delete hitForRec; - delete trackParam; - newMuonTrackRef->Delete(); - delete newMuonTrackRef; - } //_____________________________________________________________________________ -void AliMUONRecoCheck::ReconstructibleTracks() +void AliMUONRecoCheck::MakeReconstructibleTracks() { - // calculate the number of reconstructible tracks + /// Isolate the reconstructible tracks + fRecoTrackRefStore = new AliMUONTrackStoreV1(); + // create iterator on trackRef + TIter next(fTrackRefStore->CreateIterator()); + + // loop over trackRef AliMUONTrack* track; - TClonesArray* hitForRecAtHit = NULL; - AliMUONHitForRec* hitForRec; - Float_t zRec; - Int_t nTrackRef, nTrackHits; - Int_t isChamberInTrack[10]; - Int_t iChamber = 0; - Bool_t isTrackOK = kTRUE; - - fReconstructibleTracks = 0; - - nTrackRef = fMuonTrackRef->GetEntriesFast(); - for (Int_t index = 0; index < nTrackRef; index++) { - track = (AliMUONTrack*)fMuonTrackRef->At(index); - hitForRecAtHit = track->GetHitForRecAtHit(); - nTrackHits = hitForRecAtHit->GetEntriesFast(); - for (Int_t ch = 0; ch < 10; ch++) isChamberInTrack[ch] = 0; - - for ( Int_t iHit = 0; iHit < nTrackHits; iHit++) { - hitForRec = (AliMUONHitForRec*) hitForRecAtHit->At(iHit); - zRec = hitForRec->GetZ(); - iChamber = hitForRec->GetChamberNumber(); - if (iChamber < 0 || iChamber > 10) continue; - isChamberInTrack[iChamber] = 1; - - } - // track is reconstructible if the particle is crossing every tracking chambers - isTrackOK = kTRUE; - for (Int_t ch = 0; ch < 10; ch++) { - if (!isChamberInTrack[ch]) isTrackOK = kFALSE; - } if (isTrackOK) fReconstructibleTracks++; - if (!isTrackOK) fMuonTrackRef->Remove(track); // remove non reconstructible tracks + while ( ( track = static_cast(next()) ) ) { + + Bool_t* chamberInTrack = new Bool_t(AliMUONConstants::NTrackingCh()); + for (Int_t iCh = 0; iCh < AliMUONConstants::NTrackingCh(); iCh++) chamberInTrack[iCh] = kFALSE; + + // loop over trackRef's hits to get hit chambers + Int_t nTrackHits = track->GetNClusters(); + for (Int_t iHit = 0; iHit < nTrackHits; iHit++) { + AliMUONVCluster* hit = ((AliMUONTrackParam*) track->GetTrackParamAtCluster()->UncheckedAt(iHit))->GetClusterPtr(); + chamberInTrack[hit->GetChamberId()] = kTRUE; + } + + // track is reconstructible if the particle is depositing a hit + // in the following chamber combinations: + Bool_t trackOK = kTRUE; + if (!chamberInTrack[0] && !chamberInTrack[1]) trackOK = kFALSE; + if (!chamberInTrack[2] && !chamberInTrack[3]) trackOK = kFALSE; + if (!chamberInTrack[4] && !chamberInTrack[5]) trackOK = kFALSE; + Int_t nHitsInLastStations = 0; + for (Int_t iCh = 6; iCh < AliMUONConstants::NTrackingCh(); iCh++) + if (chamberInTrack[iCh]) nHitsInLastStations++; + if(nHitsInLastStations < 3) trackOK = kFALSE; + + // Add reconstructible tracks to fRecoTrackRefStore + if (trackOK) fRecoTrackRefStore->Add(*track); + + delete [] chamberInTrack; } - fMuonTrackRef->Compress(); -} +}