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