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
  58 ClassImp(AliMUONRecoEvent)
  59
  60 //-------------------------------------------------------------------
  61 AliMUONRecoEvent::AliMUONRecoEvent(Int_t eventNo)
  62   : TObject()
  63 {
  64 // Reconstructed event constructor
  65    fTracks      = new TClonesArray("AliMUONRecoTrack",200);
  66    fNevr        = eventNo;
  67    fNtracks = 0;
  68 }
  69
  70 //-------------------------------------------------------------------
  71 AliMUONRecoEvent::AliMUONRecoEvent(const AliMUONRecoEvent& rhs)
  72   : TObject(rhs)
  73 {
  74 // Protected copy constructor
  75
  76   Fatal("AliMUONRecoEventModule", "Not implemented.");
  77 }
  78
  79 //-------------------------------------------------------------------
  80 AliMUONRecoEvent::~AliMUONRecoEvent()
  81 {
  82 // Destructor of AliMUONRecoEvent
  83    fTracks->Delete();
  84    delete fTracks;
  85    fTracks = 0;
  86 }
  87
  88 //-------------------------------------------------------------------
  89 AliMUONRecoEvent&
  90 AliMUONRecoEvent::operator=(const AliMUONRecoEvent& rhs)
  91 {
  92 // Protected assignement operator
  93
  94   if (this == &rhs) return *this;
  95
  96   Fatal("operator=", "Not implemented.");
  97
  98   return *this;
  99 }
 100
 101 //-------------------------------------------------------------------
 102 AliMUONRecoTrack* AliMUONRecoEvent::AddEmptyTrack()
 103 {
 104 // Add a empty AliMUONRecoTrackObject to the track list
 105    TClonesArray &dumptracks = *fTracks;
 106    return (new(dumptracks[fNtracks++])AliMUONRecoTrack(kTRUE));
 107 }
 108
 109 //-------------------------------------------------------------------
 110 void AliMUONRecoEvent::Clear(Option_t * /*option*/)
 111 {
 112 // Clears all track pointers from the list
 113 //   fTracks->Clear(option);
 114    fTracks->Delete();
 115    fNtracks=0;
 116 }
 117
 118 //-------------------------------------------------------------------
 119 void AliMUONRecoEvent::EventInfo()
 120 {
 121 // Prints reconstructed event information
 122    cout << "*********************Reco Dumper**************************" << endl;
 123    cout << "Event number : " << fNevr << endl;
 124    cout << "   Number of tracks : " << fNtracks << endl;
 125    AliMUONRecoTrack *currentTrack =0;
 126    Int_t trackIndex = 0;
 127    for(trackIndex=0; trackIndex<fNtracks; trackIndex++) {
 128       currentTrack = (AliMUONRecoTrack*)fTracks->UncheckedAt(trackIndex);
 129       cout << "Track : " << trackIndex << endl;
 130       cout << "   Sign : " << currentTrack->GetSign() << endl;
 131       cout << "   Vertex position    : " << currentTrack->GetVertexPos() << endl;
 132       Double_t momreco[3];
 133       for (Int_t axis=0; axis<3; axis++) {
 134          momreco[axis] = currentTrack->GetMomReconstr(axis);
 135       }
 136       cout << "   Reconstructed mom. : " << "Px=" << momreco[0] << "Py=" << momreco[1] << "Pz=" << momreco[2] << endl;
 137       cout << "   Chi squared        : " << currentTrack->GetChi2r() << endl;
 138       cout << "   Hits positions     : " << endl;
 139       Double_t xhit, yhit, zhit;
 140       for (Int_t chamber=0; chamber<10; chamber++) {
 141          xhit = currentTrack->GetPosX(chamber);
 142          yhit = currentTrack->GetPosY(chamber);
 143          zhit = currentTrack->GetPosZ(chamber);
 144 //         cout <<"      chamber" << chamber << " X=" << xhit << " Y=" << yhit << " Z=" << zhit << endl;
 145       }
 146    }
 147    cout << "**********************************************************" << endl;
 148 }
 149
 150 //-------------------------------------------------------------------
 151 Bool_t AliMUONRecoEvent::MakeDumpTracks(Int_t muons, TClonesArray *tracksPtr,
 152   AliMUONEventReconstructor *EventReco)
 153 {
 154 // This method takes the pointer of the list of reconstructed tracks from
 155 // AliMUONEventReconstructor and fill the reconstructed AliMUONRecoEvent
 156 // fields.
 157
 158         cout << "Enter MakeDumpTracks..." << endl;
 159    Int_t nTracks = tracksPtr->GetEntriesFast();
 160    cout << "nTracks = "<< nTracks << endl;
 161    if (nTracks == 0) {
 162       cout << "AliMUONRecoEvent::MakeDumpTracks: Number of tracks is zero !" << endl;
 163       //AZ return kFALSE;
 164    }
 165    cout << tracksPtr << endl;
 166    if (!tracksPtr) {
 167       cout << "AliMUONRecoEvent::MakeDumpTracks() : You should call SetRecoTracksPtr() first..." << endl;
 168       return kFALSE;
 169    }
 170         // Get event number
 171    Int_t noEvent = gAlice->GetHeader()->GetEvent();
 172    cout << "noEvent = "<< nTracks << endl;
 173    tracksPtr->Compress();  // simple loop
 174    AliMUONRecoTrack *currentTrack;
 175    Int_t trackIndex, nTrackHits = 0;
 176    Double_t z, pYZ, bendingSlope, nonBendingSlope;
 177    Double_t pX, pY, pZ;                 // reconstructed momentum components
 178    Int_t isign, flag=0;         // charge sign, flag of reconstructed track
 179    Double_t alpha, beta;
 180    TObjArray *hitsOnTrack = 0;
 181    AliMUONTrackHit *trackHit = 0;
 182    AliMUONTrack *track = 0;
 183    AliMUONTrackK *trackK = 0;
 184    TMatrixD *trackParamK; //AZnon
 185    AliMUONTrackParam *trackParam = 0;
 186    // Fill event number and number of tracks
 187    fNevr = noEvent;
 188    fMuons = muons; //AZ - number of muons within acceptance
 189    // Loop over reconstructed tracks
 190    for (trackIndex=0; trackIndex<nTracks; trackIndex++) {
 191       cout << " trackIndex = " << trackIndex << endl;
 192       currentTrack = AddEmptyTrack();
 193       cout << " currentTrack = " << currentTrack << endl;
 194
 195       if (EventReco->GetTrackMethod() == 2) { // Kalman
 196
 197         trackK = (AliMUONTrackK*) ((*tracksPtr)[trackIndex]);
 198         nTrackHits = trackK->GetNTrackHits();
 199         trackParamK = trackK->GetTrackParameters();
 200         isign = Int_t(TMath::Sign(1., (*trackParamK)(4,0)));
 201         z = trackK->GetZ();
 202         alpha = (*trackParamK)(2,0);
 203         beta = (*trackParamK)(3,0);
 204         pYZ = TMath::Cos(beta)/TMath::Abs((*trackParamK)(4,0));
 205         pZ = pYZ/TMath::Cos(alpha);
 206         pX = TMath::Sin(beta)/TMath::Abs((*trackParamK)(4,0));
 207         pY = pYZ*TMath::Sin(alpha);
 208
 209         currentTrack->SetVertexPos(z);
 210         currentTrack->SetMomReconstr(pX,pY,pZ);
 211         currentTrack->SetSign(isign);
 212         currentTrack->SetChi2r(trackK->GetTrackQuality());
 213
 214         // Check hits on the track
 215         hitsOnTrack = trackK->GetHitOnTrack();
 216         Float_t signal = 0;
 217         Float_t tht = 0;
 218         for (int ihit = 0; ihit < nTrackHits; ihit++) {
 219           signal += ((AliMUONHitForRec*)((*hitsOnTrack)[ihit]))->GetGeantSignal();
 220           tht += TMath::Min (1,((AliMUONHitForRec*)((*hitsOnTrack)[ihit]))->GetTHTrack());
 221         }
 222         signal /= nTrackHits;
 223         tht /= nTrackHits;
 224         flag = 0;
 225         if (TMath::Nint(signal) > 0) { // signal muon
 226           for (int ihit = 0; ihit < nTrackHits ; ihit++) {
 227             if (((AliMUONHitForRec*)((*hitsOnTrack)[ihit]))->GetTHTrack() != TMath::Nint(tht)) flag++;
 228           }
 229         } else flag = -9; // background track
 230         //cout << TMath::Nint(signal) << " " << TMath::Nint(tht) << " " << recTrackNt->fFlag << endl;
 231         currentTrack->SetFlag(flag);
 232       } else { // default tracking
 233
 234         track = (AliMUONTrack*) ((*tracksPtr)[trackIndex]);
 235         nTrackHits = track->GetNTrackHits();
 236         // track parameters at Vertex
 237         trackParam = track->GetTrackParamAtVertex();
 238         bendingSlope = trackParam->GetBendingSlope();
 239         nonBendingSlope = trackParam->GetNonBendingSlope();
 240
 241         z = trackParam->GetZ();
 242         pYZ = 1/TMath::Abs(trackParam->GetInverseBendingMomentum());
 243         pZ = pYZ/TMath::Sqrt(1+bendingSlope*bendingSlope);
 244         pX = pZ * nonBendingSlope;
 245         pY = pZ * bendingSlope;
 246
 247         if (trackParam->GetInverseBendingMomentum()<0) isign=-1; else isign=1;
 248         currentTrack->SetVertexPos(z);
 249         currentTrack->SetMomReconstr(pX,pY,pZ);
 250         currentTrack->SetSign(isign);
 251         //         currentTrack->SetChi2r(trackParam->GetChi2());
 252         currentTrack->SetChi2r(0);
 253
 254         // Check hits on the track
 255         hitsOnTrack = track->GetTrackHitsPtr();
 256         Float_t signal = 0;
 257         Float_t tht = 0;
 258         AliMUONHitForRec *hitForRec = 0;
 259         for (int ihit = 0; ihit < nTrackHits; ihit++) {
 260           hitForRec = ((AliMUONTrackHit*)(*hitsOnTrack)[ihit])->GetHitForRecPtr();
 261           signal += hitForRec->GetGeantSignal();
 262           tht += TMath::Min (1,hitForRec->GetTHTrack());
 263         }
 264         signal /= nTrackHits;
 265         tht /= nTrackHits;
 266         flag = 0;
 267         if (TMath::Nint(signal) > 0) { // signal muon
 268           for (int ihit = 0; ihit < nTrackHits ; ihit++) {
 269             hitForRec = ((AliMUONTrackHit*)(*hitsOnTrack)[ihit])->GetHitForRecPtr();
 270             if (hitForRec->GetTHTrack() != TMath::Nint(tht)) flag++;
 271           }
 272         } else flag = -9; // background track
 273         //cout << TMath::Nint(signal) << " " << TMath::Nint(tht) << " " << recTrackNt->fFlag << endl;
 274         currentTrack->SetFlag(flag);
 275       }
 276
 277       Double_t xhit,yhit,zhit;
 278       // Loop over track hits
 279       for (Int_t trackHitIndex = 0; trackHitIndex < nTrackHits; trackHitIndex++) {
 280         if (EventReco->GetTrackMethod() == 2) { // Kalman
 281           xhit = ((AliMUONHitForRec*)((*hitsOnTrack)[trackHitIndex]))->GetNonBendingCoor();
 282           yhit = ((AliMUONHitForRec*)((*hitsOnTrack)[trackHitIndex]))->GetBendingCoor();
 283           zhit = ((AliMUONHitForRec*)((*hitsOnTrack)[trackHitIndex]))->GetZ();
 284         } else {
 285           trackHit = (AliMUONTrackHit*) (*(track->GetTrackHitsPtr()))[trackHitIndex];
 286           xhit = trackHit->GetHitForRecPtr()->GetNonBendingCoor();
 287           yhit = trackHit->GetHitForRecPtr()->GetBendingCoor();
 288           zhit = trackHit->GetHitForRecPtr()->GetZ();
 289         }
 290         if (trackHitIndex >= 0 && trackHitIndex < 10) {
 291           currentTrack->SetHitPosition(trackHitIndex,xhit,yhit,zhit);
 292         } else { cout << "track " << trackIndex << " hit out of range" << endl;}
 293       }
 294
 295    }
 296    cout << "Leave MakeDumpTracks..." << endl;
 297    return kTRUE;
 298 }
 299
 300 //-------------------------------------------------------------------
 301 void AliMUONRecoEvent::Streamer(TBuffer &R__b)
 302 {
 303 // Streams an object of class AliMUONRecoEvent
 304    if (R__b.IsReading()) {
 305       fTracks->Clear();
 306       AliMUONRecoEvent::Class()->ReadBuffer(R__b, this);
 307    } else {
 308       cout << "...writing event to file...\n";
 309       AliMUONRecoEvent::Class()->WriteBuffer(R__b, this);
 310    }
 311 }