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 /*
  17 $Log$
  18
  19 Revision 1.7  2002/10/23 07:24:56  alibrary
  20 Introducing Riostream.h
  21
  22 Revision 1.6  2001/05/16 14:57:17  alibrary
  23 New files for folders and Stack
  24
  25 Revision 1.5  2001/05/11 07:59:03  hristov
  26 Index corrected
  27
  28 Revision 1.4  2001/01/26 20:00:53  hristov
  29 Major upgrade of AliRoot code
  30
  31 Revision 1.3  2000/12/21 17:51:54  morsch
  32 RN3 violations corrected
  33
  34 Revision 1.2  2000/11/23 10:09:38  gosset
  35 Bug correction in AliMUONRecoDisplay.
  36 Copyright, $Log$
  37 Copyright, Revision 1.7  2002/10/23 07:24:56  alibrary
  38 Copyright, Introducing Riostream.h
  39 Copyright,
  40 Copyright, Revision 1.6  2001/05/16 14:57:17  alibrary
  41 Copyright, New files for folders and Stack
  42 Copyright,
  43 Copyright, Revision 1.5  2001/05/11 07:59:03  hristov
  44 Copyright, Index corrected
  45 Copyright,
  46 Copyright, Revision 1.4  2001/01/26 20:00:53  hristov
  47 Copyright, Major upgrade of AliRoot code
  48 Copyright,
  49 Copyright, Revision 1.3  2000/12/21 17:51:54  morsch
  50 Copyright, RN3 violations corrected
  51 Copyright,, $Id$, comments at the right place for automatic documentation,
  52 in AliMUONRecoEvent and AliMUONRecoDisplay
  53
  54
  55 */
  56
  57 //Authors: Mihaela Gheata, Andrei Gheata 09/10/00
  58 ////////////////////////////////////////////////////////////////////
  59 //                                                                //
  60 // AliMUONRecoEvent, AliMUONRecoTrack (and AliMUONRecoDisplay)    //
  61 //                                                                //
  62 // Theses classes are used to store and retrieve                  //
  63 // MUON reconstructed events.                                     //
  64 // The corresponding tree is constructed and filled               //
  65 // during the FillEvent() method of AliMUONEventReconstructor,    //
  66 // when all reconstruction information is available.              //
  67 //                                                                //
  68 ////////////////////////////////////////////////////////////////////
  69 ////////////////////////////////////////////////////////////////////
  70 //                                                                //
  71 // AliMUONRecoEvent                                               //
  72 //                                                                //
  73 // This class handles an array of reconstructed tracks.           //
  74 // It provides :                                                  //
  75 //      - filling the tracks array according to the information   //
  76 //        stored in AliMUONEventReconstructor class ;             //
  77 //      - printing event and track informations : event numer,    //
  78 //        number of tracks, hits positions, reconstr. mometum.    //
  79 //                                                                //
  80 ////////////////////////////////////////////////////////////////////
  81 ////////////////////////////////////////////////////////////////////
  82 //                                                                //
  83 // AliMUONRecoTrack                                               //
  84 //                                                                //
  85 // This class represents a reconstructed muon track               //
  86 // in the tree of reconstructed events.                           //
  87 //                                                                //
  88 ////////////////////////////////////////////////////////////////////
  89
  90 #include <Riostream.h>
  91 #include <AliRun.h>
  92 #include <TClonesArray.h>
  93 #include <TClass.h>
  94
  95 #include <TFile.h>
  96 #include <TMatrixD.h>
  97 #include <TParticle.h>
  98
  99 #include "AliMUONRecoEvent.h"
 100 #include "AliMUONEventReconstructor.h"
 101 #include "AliMUONTrack.h"
 102 #include "AliMUONTrackK.h"
 103 #include "AliMUONTrackParam.h"
 104 #include "AliMUONHitForRec.h"
 105 #include "AliMUONTrackHit.h"
 106 #include "AliHeader.h"
 107
 108 ClassImp(AliMUONRecoTrack)
 109 ClassImp(AliMUONRecoEvent)
 110
 111 //-------------------------------------------------------------------
 112 AliMUONRecoEvent::AliMUONRecoEvent(Int_t eventNo)
 113 {
 114 // Reconstructed event constructor
 115    fTracks      = new TClonesArray("AliMUONRecoTrack",200);
 116    fNevr        = eventNo;
 117    fNtracks = 0;
 118 }
 119
 120 //-------------------------------------------------------------------
 121 AliMUONRecoEvent::~AliMUONRecoEvent()
 122 {
 123 // Destructor of AliMUONRecoEvent
 124    fTracks->Delete();
 125    delete fTracks;
 126    fTracks = 0;
 127 }
 128
 129 //-------------------------------------------------------------------
 130 AliMUONRecoTrack* AliMUONRecoEvent::AddEmptyTrack()
 131 {
 132 // Add a empty AliMUONRecoTrackObject to the track list
 133    TClonesArray &dumptracks = *fTracks;
 134    return (new(dumptracks[fNtracks++])AliMUONRecoTrack(kTRUE));
 135 }
 136
 137 //-------------------------------------------------------------------
 138 void AliMUONRecoEvent::Clear(Option_t *option)
 139 {
 140 // Clears all track pointers from the list
 141 //   fTracks->Clear(option);
 142    fTracks->Delete();
 143    fNtracks=0;
 144 }
 145
 146 //-------------------------------------------------------------------
 147 void AliMUONRecoEvent::EventInfo()
 148 {
 149 // Prints reconstructed event information
 150    cout << "*********************Reco Dumper**************************" << endl;
 151    cout << "Event number : " << fNevr << endl;
 152    cout << "   Number of tracks : " << fNtracks << endl;
 153    AliMUONRecoTrack *currentTrack =0;
 154    Int_t trackIndex = 0;
 155    for(trackIndex=0; trackIndex<fNtracks; trackIndex++) {
 156       currentTrack = (AliMUONRecoTrack*)fTracks->UncheckedAt(trackIndex);
 157       cout << "Track : " << trackIndex << endl;
 158       cout << "   Sign : " << currentTrack->GetSign() << endl;
 159       cout << "   Vertex position    : " << currentTrack->GetVertexPos() << endl;
 160       Double_t momreco[3];
 161       for (Int_t axis=0; axis<3; axis++) {
 162          momreco[axis] = currentTrack->GetMomReconstr(axis);
 163       }
 164       cout << "   Reconstructed mom. : " << "Px=" << momreco[0] << "Py=" << momreco[1] << "Pz=" << momreco[2] << endl;
 165       cout << "   Chi squared        : " << currentTrack->GetChi2r() << endl;
 166       cout << "   Hits positions     : " << endl;
 167       Double_t xhit, yhit, zhit;
 168       for (Int_t chamber=0; chamber<10; chamber++) {
 169          xhit = currentTrack->GetPosX(chamber);
 170          yhit = currentTrack->GetPosY(chamber);
 171          zhit = currentTrack->GetPosZ(chamber);
 172 //         cout <<"      chamber" << chamber << " X=" << xhit << " Y=" << yhit << " Z=" << zhit << endl;
 173       }
 174    }
 175    cout << "**********************************************************" << endl;
 176 }
 177
 178 //-------------------------------------------------------------------
 179 Bool_t AliMUONRecoEvent::MakeDumpTracks(Int_t muons, TClonesArray *tracksPtr,
 180   AliMUONEventReconstructor *EventReco)
 181 {
 182 // This method takes the pointer of the list of reconstructed tracks from
 183 // AliMUONEventReconstructor and fill the reconstructed AliMUONRecoEvent
 184 // fields.
 185
 186         cout << "Enter MakeDumpTracks..." << endl;
 187    Int_t nTracks = tracksPtr->GetEntriesFast();
 188    cout << "nTracks = "<< nTracks << endl;
 189    if (nTracks == 0) {
 190       cout << "AliMUONRecoEvent::MakeDumpTracks: Number of tracks is zero !" << endl;
 191       //AZ return kFALSE;
 192    }
 193    cout << tracksPtr << endl;
 194    if (!tracksPtr) {
 195       cout << "AliMUONRecoEvent::MakeDumpTracks() : You should call SetRecoTracksPtr() first..." << endl;
 196       return kFALSE;
 197    }
 198         // Get event number
 199    Int_t noEvent = gAlice->GetHeader()->GetEvent();
 200    cout << "noEvent = "<< nTracks << endl;
 201    tracksPtr->Compress();  // simple loop
 202    AliMUONRecoTrack *currentTrack;
 203    Int_t trackIndex, nTrackHits = 0;
 204    Double_t z, pYZ, bendingSlope, nonBendingSlope;
 205    Double_t pX, pY, pZ;                 // reconstructed momentum components
 206    Int_t isign, flag=0;         // charge sign, flag of reconstructed track
 207    Double_t alpha, beta;
 208    TObjArray *hitsOnTrack = 0;
 209    AliMUONTrackHit *trackHit = 0;
 210    AliMUONTrack *track = 0;
 211    AliMUONTrackK *trackK = 0;
 212    TMatrixD *trackParamK; //AZnon
 213    AliMUONTrackParam *trackParam = 0;
 214    // Fill event number and number of tracks
 215    fNevr = noEvent;
 216    fMuons = muons; //AZ - number of muons within acceptance
 217    // Loop over reconstructed tracks
 218    for (trackIndex=0; trackIndex<nTracks; trackIndex++) {
 219       cout << " trackIndex = " << trackIndex << endl;
 220       currentTrack = AddEmptyTrack();
 221       cout << " currentTrack = " << currentTrack << endl;
 222
 223       if (EventReco->GetTrackMethod() == 2) { // Kalman
 224
 225         trackK = (AliMUONTrackK*) ((*tracksPtr)[trackIndex]);
 226         nTrackHits = trackK->GetNTrackHits();
 227         trackParamK = trackK->GetTrackParameters();
 228         isign = Int_t(TMath::Sign(1., (*trackParamK)(4,0)));
 229         z = trackK->GetZ();
 230         alpha = (*trackParamK)(2,0);
 231         beta = (*trackParamK)(3,0);
 232         pYZ = TMath::Cos(beta)/TMath::Abs((*trackParamK)(4,0));
 233         pZ = pYZ/TMath::Cos(alpha);
 234         pX = TMath::Sin(beta)/TMath::Abs((*trackParamK)(4,0));
 235         pY = pYZ*TMath::Sin(alpha);
 236
 237         currentTrack->SetVertexPos(z);
 238         currentTrack->SetMomReconstr(pX,pY,pZ);
 239         currentTrack->SetSign(isign);
 240         currentTrack->SetChi2r(trackK->GetTrackQuality());
 241
 242         // Check hits on the track
 243         hitsOnTrack = trackK->GetHitOnTrack();
 244         Float_t signal = 0;
 245         Float_t tht = 0;
 246         for (int ihit = 0; ihit < nTrackHits; ihit++) {
 247           signal += ((AliMUONHitForRec*)((*hitsOnTrack)[ihit]))->GetGeantSignal();
 248           tht += TMath::Min (1,((AliMUONHitForRec*)((*hitsOnTrack)[ihit]))->GetTHTrack());
 249         }
 250         signal /= nTrackHits;
 251         tht /= nTrackHits;
 252         flag = 0;
 253         if (TMath::Nint(signal) > 0) { // signal muon
 254           for (int ihit = 0; ihit < nTrackHits ; ihit++) {
 255             if (((AliMUONHitForRec*)((*hitsOnTrack)[ihit]))->GetTHTrack() != TMath::Nint(tht)) flag++;
 256           }
 257         } else flag = -9; // background track
 258         //cout << TMath::Nint(signal) << " " << TMath::Nint(tht) << " " << recTrackNt->fFlag << endl;
 259         currentTrack->SetFlag(flag);
 260       } else { // default tracking
 261
 262         track = (AliMUONTrack*) ((*tracksPtr)[trackIndex]);
 263         nTrackHits = track->GetNTrackHits();
 264         // track parameters at Vertex
 265         trackParam = track->GetTrackParamAtVertex();
 266         bendingSlope = trackParam->GetBendingSlope();
 267         nonBendingSlope = trackParam->GetNonBendingSlope();
 268
 269         z = trackParam->GetZ();
 270         pYZ = 1/TMath::Abs(trackParam->GetInverseBendingMomentum());
 271         pZ = pYZ/TMath::Sqrt(1+bendingSlope*bendingSlope);
 272         pX = pZ * nonBendingSlope;
 273         pY = pZ * bendingSlope;
 274
 275         if (trackParam->GetInverseBendingMomentum()<0) isign=-1; else isign=1;
 276         currentTrack->SetVertexPos(z);
 277         currentTrack->SetMomReconstr(pX,pY,pZ);
 278         currentTrack->SetSign(isign);
 279         //         currentTrack->SetChi2r(trackParam->GetChi2());
 280         currentTrack->SetChi2r(0);
 281
 282         // Check hits on the track
 283         hitsOnTrack = track->GetTrackHitsPtr();
 284         Float_t signal = 0;
 285         Float_t tht = 0;
 286         AliMUONHitForRec *hitForRec = 0;
 287         for (int ihit = 0; ihit < nTrackHits; ihit++) {
 288           hitForRec = ((AliMUONTrackHit*)(*hitsOnTrack)[ihit])->GetHitForRecPtr();
 289           signal += hitForRec->GetGeantSignal();
 290           tht += TMath::Min (1,hitForRec->GetTHTrack());
 291         }
 292         signal /= nTrackHits;
 293         tht /= nTrackHits;
 294         flag = 0;
 295         if (TMath::Nint(signal) > 0) { // signal muon
 296           for (int ihit = 0; ihit < nTrackHits ; ihit++) {
 297             hitForRec = ((AliMUONTrackHit*)(*hitsOnTrack)[ihit])->GetHitForRecPtr();
 298             if (hitForRec->GetTHTrack() != TMath::Nint(tht)) flag++;
 299           }
 300         } else flag = -9; // background track
 301         //cout << TMath::Nint(signal) << " " << TMath::Nint(tht) << " " << recTrackNt->fFlag << endl;
 302         currentTrack->SetFlag(flag);
 303       }
 304
 305       Double_t xhit,yhit,zhit;
 306       // Loop over track hits
 307       for (Int_t trackHitIndex = 0; trackHitIndex < nTrackHits; trackHitIndex++) {
 308         if (EventReco->GetTrackMethod() == 2) { // Kalman
 309           xhit = ((AliMUONHitForRec*)((*hitsOnTrack)[trackHitIndex]))->GetNonBendingCoor();
 310           yhit = ((AliMUONHitForRec*)((*hitsOnTrack)[trackHitIndex]))->GetBendingCoor();
 311           zhit = ((AliMUONHitForRec*)((*hitsOnTrack)[trackHitIndex]))->GetZ();
 312         } else {
 313           trackHit = (AliMUONTrackHit*) (*(track->GetTrackHitsPtr()))[trackHitIndex];
 314           xhit = trackHit->GetHitForRecPtr()->GetNonBendingCoor();
 315           yhit = trackHit->GetHitForRecPtr()->GetBendingCoor();
 316           zhit = trackHit->GetHitForRecPtr()->GetZ();
 317         }
 318         if (trackHitIndex >= 0 && trackHitIndex < 10) {
 319           currentTrack->SetHitPosition(trackHitIndex,xhit,yhit,zhit);
 320         } else { cout << "track " << trackIndex << " hit out of range" << endl;}
 321       }
 322
 323    }
 324    cout << "Leave MakeDumpTracks..." << endl;
 325    return kTRUE;
 326 }
 327
 328 //-------------------------------------------------------------------
 329 void AliMUONRecoEvent::Streamer(TBuffer &R__b)
 330 {
 331 // Streams an object of class AliMUONRecoEvent
 332    if (R__b.IsReading()) {
 333       fTracks->Clear();
 334       AliMUONRecoEvent::Class()->ReadBuffer(R__b, this);
 335    } else {
 336       cout << "...writing event to file...\n";
 337       AliMUONRecoEvent::Class()->WriteBuffer(R__b, this);
 338    }
 339 }
 340
 341 //-------------------------------------------------------------------
 342 AliMUONRecoTrack::AliMUONRecoTrack(Bool_t active)
 343 {
 344 //Constructor of AliMUONRecoTrack
 345    fSign  = 0;
 346    fZvr   = 0.0;
 347    fChi2r = 0.0;
 348    if (active) {
 349         for (Int_t axis=0; axis<3; axis++) {
 350         fPr[axis] = 0.0;
 351         }
 352         for (Int_t chamber=0; chamber<10; chamber++) {
 353         fPosX[chamber] = 0.0;
 354         fPosY[chamber] = 0.0;
 355         fPosZ[chamber] = 0.0;
 356         }
 357    }
 358 }
 359
 360 //-------------------------------------------------------------------
 361 const Double_t AliMUONRecoTrack::Phi()
 362 {
 363 // Return trach phi angle
 364         return TMath::ATan2(fPr[2], fPr[1]);
 365 }
 366
 367 //-------------------------------------------------------------------
 368 const Double_t AliMUONRecoTrack::Theta()
 369 {
 370 // Return trach theta angle
 371    return TMath::ACos(fPr[2] / P());
 372 }
 373
 374 //-------------------------------------------------------------------
 375 void AliMUONRecoTrack::SetMomReconstr(Double_t px, Double_t py, Double_t pz)
 376 {
 377 // Set the track reconstructed momentum
 378    fPr[0] = px;
 379    fPr[1] = py;
 380    fPr[2] = pz;
 381 }
 382
 383 //-------------------------------------------------------------------
 384 void AliMUONRecoTrack::SetHitPosition(Int_t chamber, Double_t x, Double_t y, Double_t z)
 385 {
 386 // Set hit coordinates in chamber[0..9]
 387    fPosX[chamber] = x;
 388    fPosY[chamber] = y;
 389    fPosZ[chamber] = z;
 390 }
 391 //-------------------------------------------------------------------
 392 void AliMUONRecoTrack::TrackInfo()
 393 {
 394 // Prints momentum info for this track
 395    cout << "Px=" << GetMomReconstr(0) << " Py=" << GetMomReconstr(1) <<
 396           " Pz=" << GetMomReconstr(2) << " P=" << P() << endl;
 397 }