MUON/AliMUONRecoEvent.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /* $Id$ */
  17
  18 //Authors: Mihaela Gheata, Andrei Gheata 09/10/00
  19 ////////////////////////////////////////////////////////////////////
  20 //                                                                //
  21 // AliMUONRecoEvent, AliMUONRecoTrack (and AliMUONRecoDisplay)    //
  22 //                                                                //
  23 // Theses classes are used to store and retrieve                  //
  24 // MUON reconstructed events.                                     //
  25 // The corresponding tree is constructed and filled               //
  26 // during the FillEvent() method of AliMUONEventReconstructor,    //
  27 // when all reconstruction information is available.              //
  28 //                                                                //
  29 ////////////////////////////////////////////////////////////////////
  30 ////////////////////////////////////////////////////////////////////
  31 //                                                                //
  32 // AliMUONRecoEvent                                               //
  33 //                                                                //
  34 // This class handles an array of reconstructed tracks.           //
  35 // It provides :                                                  //
  36 //      - filling the tracks array according to the information   //
  37 //        stored in AliMUONEventReconstructor class ;             //
  38 //      - printing event and track informations : event numer,    //
  39 //        number of tracks, hits positions, reconstr. mometum.    //
  40 //                                                                //
  41 ////////////////////////////////////////////////////////////////////
  42
  43 #include <Riostream.h>
  44 #include <TClonesArray.h>
  45 #include <TMatrixD.h>
  46
  47 #include "AliMUONRecoEvent.h"
  48 #include "AliMUONRecoTrack.h"
  49 #include "AliMUONEventReconstructor.h"
  50 #include "AliMUONTrack.h"
  51 #include "AliMUONTrackK.h"
  52 #include "AliMUONTrackParam.h"
  53 #include "AliMUONHitForRec.h"
  54 #include "AliMUONTrackHit.h"
  55 #include "AliRun.h"
  56 #include "AliHeader.h"
  57 #include "AliLog.h"
  58
  59 ClassImp(AliMUONRecoEvent)
  60
  61 //-------------------------------------------------------------------
  62 AliMUONRecoEvent::AliMUONRecoEvent(Int_t eventNo)
  63   : TObject()
  64 {
  65 // Reconstructed event constructor
  66    fTracks      = new TClonesArray("AliMUONRecoTrack",200);
  67    fNevr        = eventNo;
  68    fNtracks = 0;
  69 }
  70
  71 //-------------------------------------------------------------------
  72 AliMUONRecoEvent::AliMUONRecoEvent(const AliMUONRecoEvent& rhs)
  73   : TObject(rhs)
  74 {
  75 // Protected copy constructor
  76
  77   AliFatal("Not implemented.");
  78 }
  79
  80 //-------------------------------------------------------------------
  81 AliMUONRecoEvent::~AliMUONRecoEvent()
  82 {
  83 // Destructor of AliMUONRecoEvent
  84    fTracks->Delete();
  85    delete fTracks;
  86    fTracks = 0;
  87 }
  88
  89 //-------------------------------------------------------------------
  90 AliMUONRecoEvent&
  91 AliMUONRecoEvent::operator=(const AliMUONRecoEvent& rhs)
  92 {
  93 // Protected assignement operator
  94
  95   if (this == &rhs) return *this;
  96
  97   AliFatal("Not implemented.");
  98
  99   return *this;
 100 }
 101
 102 //-------------------------------------------------------------------
 103 AliMUONRecoTrack* AliMUONRecoEvent::AddEmptyTrack()
 104 {
 105 // Add a empty AliMUONRecoTrackObject to the track list
 106    TClonesArray &dumptracks = *fTracks;
 107    return (new(dumptracks[fNtracks++])AliMUONRecoTrack(kTRUE));
 108 }
 109
 110 //-------------------------------------------------------------------
 111 void AliMUONRecoEvent::Clear(Option_t * /*option*/)
 112 {
 113 // Clears all track pointers from the list
 114 //   fTracks->Clear(option);
 115    fTracks->Delete();
 116    fNtracks=0;
 117 }
 118
 119 //-------------------------------------------------------------------
 120 void AliMUONRecoEvent::EventInfo()
 121 {
 122 // Prints reconstructed event information
 123    cout << "*********************Reco Dumper**************************" << endl;
 124    cout << "Event number : " << fNevr << endl;
 125    cout << "   Number of tracks : " << fNtracks << endl;
 126    AliMUONRecoTrack *currentTrack =0;
 127    Int_t trackIndex = 0;
 128    for(trackIndex=0; trackIndex<fNtracks; trackIndex++) {
 129       currentTrack = (AliMUONRecoTrack*)fTracks->UncheckedAt(trackIndex);
 130       cout << "Track : " << trackIndex << endl;
 131       cout << "   Sign : " << currentTrack->GetSign() << endl;
 132       cout << "   Vertex position    : " << currentTrack->GetVertexPos() << endl;
 133       Double_t momreco[3];
 134       for (Int_t axis=0; axis<3; axis++) {
 135          momreco[axis] = currentTrack->GetMomReconstr(axis);
 136       }
 137       cout << "   Reconstructed mom. : " << "Px=" << momreco[0] << "Py=" << momreco[1] << "Pz=" << momreco[2] << endl;
 138       cout << "   Chi squared        : " << currentTrack->GetChi2r() << endl;
 139       cout << "   Hits positions     : " << endl;
 140       Double_t xhit, yhit, zhit;
 141       for (Int_t chamber=0; chamber<10; chamber++) {
 142          xhit = currentTrack->GetPosX(chamber);
 143          yhit = currentTrack->GetPosY(chamber);
 144          zhit = currentTrack->GetPosZ(chamber);
 145 //         cout <<"      chamber" << chamber << " X=" << xhit << " Y=" << yhit << " Z=" << zhit << endl;
 146       }
 147    }
 148    cout << "**********************************************************" << endl;
 149 }
 150
 151 //-------------------------------------------------------------------
 152 Bool_t AliMUONRecoEvent::MakeDumpTracks(Int_t muons, TClonesArray *tracksPtr,
 153   AliMUONEventReconstructor *EventReco)
 154 {
 155 // This method takes the pointer of the list of reconstructed tracks from
 156 // AliMUONEventReconstructor and fill the reconstructed AliMUONRecoEvent
 157 // fields.
 158
 159         cout << "Enter MakeDumpTracks..." << endl;
 160    Int_t nTracks = tracksPtr->GetEntriesFast();
 161    cout << "nTracks = "<< nTracks << endl;
 162    if (nTracks == 0) {
 163       cout << "AliMUONRecoEvent::MakeDumpTracks: Number of tracks is zero !" << endl;
 164       //AZ return kFALSE;
 165    }
 166    cout << tracksPtr << endl;
 167    if (!tracksPtr) {
 168       cout << "AliMUONRecoEvent::MakeDumpTracks() : You should call SetRecoTracksPtr() first..." << endl;
 169       return kFALSE;
 170    }
 171         // Get event number
 172    Int_t noEvent = gAlice->GetHeader()->GetEvent();
 173    cout << "noEvent = "<< nTracks << endl;
 174    tracksPtr->Compress();  // simple loop
 175    AliMUONRecoTrack *currentTrack;
 176    Int_t trackIndex, nTrackHits = 0;
 177    Double_t z, pYZ, bendingSlope, nonBendingSlope;
 178    Double_t pX, pY, pZ;                 // reconstructed momentum components
 179    Int_t isign, flag=0;         // charge sign, flag of reconstructed track
 180    Double_t alpha, beta;
 181    TObjArray *hitsOnTrack = 0;
 182    AliMUONTrackHit *trackHit = 0;
 183    AliMUONTrack *track = 0;
 184    AliMUONTrackK *trackK = 0;
 185    TMatrixD *trackParamK; //AZnon
 186    AliMUONTrackParam *trackParam = 0;
 187    // Fill event number and number of tracks
 188    fNevr = noEvent;
 189    fMuons = muons; //AZ - number of muons within acceptance
 190    // Loop over reconstructed tracks
 191    for (trackIndex=0; trackIndex<nTracks; trackIndex++) {
 192       cout << " trackIndex = " << trackIndex << endl;
 193       currentTrack = AddEmptyTrack();
 194       cout << " currentTrack = " << currentTrack << endl;
 195
 196       if (EventReco->GetTrackMethod() == 2) { // Kalman
 197
 198         trackK = (AliMUONTrackK*) ((*tracksPtr)[trackIndex]);
 199         nTrackHits = trackK->GetNTrackHits();
 200         trackParamK = trackK->GetTrackParameters();
 201         isign = Int_t(TMath::Sign(1., (*trackParamK)(4,0)));
 202         z = trackK->GetZ();
 203         alpha = (*trackParamK)(2,0);
 204         beta = (*trackParamK)(3,0);
 205         pYZ = TMath::Cos(beta)/TMath::Abs((*trackParamK)(4,0));
 206         pZ = pYZ/TMath::Cos(alpha);
 207         pX = TMath::Sin(beta)/TMath::Abs((*trackParamK)(4,0));
 208         pY = pYZ*TMath::Sin(alpha);
 209
 210         currentTrack->SetVertexPos(z);
 211         currentTrack->SetMomReconstr(pX,pY,pZ);
 212         currentTrack->SetSign(isign);
 213         currentTrack->SetChi2r(trackK->GetTrackQuality());
 214
 215         // Check hits on the track
 216         hitsOnTrack = trackK->GetHitOnTrack();
 217         Float_t signal = 0;
 218         Float_t tht = 0;
 219         for (int ihit = 0; ihit < nTrackHits; ihit++) {
 220           signal += ((AliMUONHitForRec*)((*hitsOnTrack)[ihit]))->GetTrackRefSignal();
 221           tht += TMath::Min (1,((AliMUONHitForRec*)((*hitsOnTrack)[ihit]))->GetTTRTrack());
 222         }
 223         signal /= nTrackHits;
 224         tht /= nTrackHits;
 225         flag = 0;
 226         if (TMath::Nint(signal) > 0) { // signal muon
 227           for (int ihit = 0; ihit < nTrackHits ; ihit++) {
 228             if (((AliMUONHitForRec*)((*hitsOnTrack)[ihit]))->GetTTRTrack() != TMath::Nint(tht)) flag++;
 229           }
 230         } else flag = -9; // background track
 231         //cout << TMath::Nint(signal) << " " << TMath::Nint(tht) << " " << recTrackNt->fFlag << endl;
 232         currentTrack->SetFlag(flag);
 233       } else { // default tracking
 234
 235         track = (AliMUONTrack*) ((*tracksPtr)[trackIndex]);
 236         nTrackHits = track->GetNTrackHits();
 237         // track parameters at Vertex
 238         trackParam = track->GetTrackParamAtVertex();
 239         bendingSlope = trackParam->GetBendingSlope();
 240         nonBendingSlope = trackParam->GetNonBendingSlope();
 241
 242         z = trackParam->GetZ();
 243         pYZ = 1/TMath::Abs(trackParam->GetInverseBendingMomentum());
 244         pZ = pYZ/TMath::Sqrt(1+bendingSlope*bendingSlope);
 245         pX = pZ * nonBendingSlope;
 246         pY = pZ * bendingSlope;
 247
 248         if (trackParam->GetInverseBendingMomentum()<0) isign=-1; else isign=1;
 249         currentTrack->SetVertexPos(z);
 250         currentTrack->SetMomReconstr(pX,pY,pZ);
 251         currentTrack->SetSign(isign);
 252         //         currentTrack->SetChi2r(trackParam->GetChi2());
 253         currentTrack->SetChi2r(0);
 254
 255         // Check hits on the track
 256         hitsOnTrack = track->GetTrackHitsPtr();
 257         Float_t signal = 0;
 258         Float_t tht = 0;
 259         AliMUONHitForRec *hitForRec = 0;
 260         for (int ihit = 0; ihit < nTrackHits; ihit++) {
 261           hitForRec = ((AliMUONTrackHit*)(*hitsOnTrack)[ihit])->GetHitForRecPtr();
 262           signal += hitForRec->GetTrackRefSignal();
 263           tht += TMath::Min (1,hitForRec->GetTTRTrack());
 264         }
 265         signal /= nTrackHits;
 266         tht /= nTrackHits;
 267         flag = 0;
 268         if (TMath::Nint(signal) > 0) { // signal muon
 269           for (int ihit = 0; ihit < nTrackHits ; ihit++) {
 270             hitForRec = ((AliMUONTrackHit*)(*hitsOnTrack)[ihit])->GetHitForRecPtr();
 271             if (hitForRec->GetTTRTrack() != TMath::Nint(tht)) flag++;
 272           }
 273         } else flag = -9; // background track
 274         //cout << TMath::Nint(signal) << " " << TMath::Nint(tht) << " " << recTrackNt->fFlag << endl;
 275         currentTrack->SetFlag(flag);
 276       }
 277
 278       Double_t xhit,yhit,zhit;
 279       // Loop over track hits
 280       for (Int_t trackHitIndex = 0; trackHitIndex < nTrackHits; trackHitIndex++) {
 281         if (EventReco->GetTrackMethod() == 2) { // Kalman
 282           xhit = ((AliMUONHitForRec*)((*hitsOnTrack)[trackHitIndex]))->GetNonBendingCoor();
 283           yhit = ((AliMUONHitForRec*)((*hitsOnTrack)[trackHitIndex]))->GetBendingCoor();
 284           zhit = ((AliMUONHitForRec*)((*hitsOnTrack)[trackHitIndex]))->GetZ();
 285         } else {
 286           trackHit = (AliMUONTrackHit*) (*(track->GetTrackHitsPtr()))[trackHitIndex];
 287           xhit = trackHit->GetHitForRecPtr()->GetNonBendingCoor();
 288           yhit = trackHit->GetHitForRecPtr()->GetBendingCoor();
 289           zhit = trackHit->GetHitForRecPtr()->GetZ();
 290         }
 291         if (trackHitIndex >= 0 && trackHitIndex < 10) {
 292           currentTrack->SetHitPosition(trackHitIndex,xhit,yhit,zhit);
 293         } else { cout << "track " << trackIndex << " hit out of range" << endl;}
 294       }
 295
 296    }
 297    cout << "Leave MakeDumpTracks..." << endl;
 298    return kTRUE;
 299 }
 300
 301 //-------------------------------------------------------------------
 302 void AliMUONRecoEvent::Streamer(TBuffer &R__b)
 303 {
 304 // Streams an object of class AliMUONRecoEvent
 305    if (R__b.IsReading()) {
 306       fTracks->Clear();
 307       AliMUONRecoEvent::Class()->ReadBuffer(R__b, this);
 308    } else {
 309       cout << "...writing event to file...\n";
 310       AliMUONRecoEvent::Class()->WriteBuffer(R__b, this);
 311    }
 312 }