MUON/AliMUON.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
  19 ///////////////////////////////////////////////
  20 //  Manager and hits classes for set:MUON     //
  21 ////////////////////////////////////////////////
  22
  23 #include "Riostream.h"
  24
  25 #include <AliPDG.h>
  26 #include <TBRIK.h>
  27 #include <TCanvas.h>
  28 #include <TDirectory.h>
  29 #include <TFile.h>
  30 #include <TGeometry.h>
  31 #include <TMinuit.h>
  32 #include <TNode.h>
  33 #include <TNtuple.h>
  34 #include <TObjArray.h>
  35 #include <TObject.h>
  36 #include <TObjectTable.h>
  37 #include <TPad.h>
  38 #include <TParticle.h>
  39 #include <TROOT.h>
  40 #include <TRandom.h>
  41 #include <TRotMatrix.h>
  42 #include <TTUBE.h>
  43 #include <TTUBE.h>
  44 #include <TTree.h>
  45 #include <TVector.h>
  46 #include <TVirtualMC.h>
  47
  48 #include "AliConst.h"
  49 #include "AliHeader.h"
  50 #include "AliHitMap.h"
  51 #include "AliLoader.h"
  52 #include "AliESD.h"
  53 #include "AliESDMuonTrack.h"
  54 #include "AliMC.h"
  55 #include "AliMUONLoader.h"
  56 #include "AliMUON.h"
  57 #include "AliMUONTriggerTrack.h"
  58 #include "AliMUONEventReconstructor.h"
  59 #include "AliMUONClusterReconstructor.h"
  60 #include "AliMUONTrack.h"
  61 #include "AliMUONTrackParam.h"
  62 #include "AliMUONChamberTrigger.h"
  63 #include "AliMUONClusterFinderAZ.h"
  64 #include "AliMUONClusterInput.h"
  65 #include "AliMUONConstants.h"
  66 #include "AliMUONDigit.h"
  67 #include "AliMUONGlobalTrigger.h"
  68 #include "AliMUONHit.h"
  69 #include "AliMUONHitMapA1.h"
  70 #include "AliMUONLocalTrigger.h"
  71 #include "AliMUONMerger.h"
  72 #include "AliMUONPadHit.h"
  73 #include "AliMUONRawCluster.h"
  74 #include "AliMUONTransientDigit.h"
  75 #include "AliMUONTriggerCircuit.h"
  76 #include "AliMUONTriggerDecision.h"
  77 #include "AliMUONVGeometryBuilder.h"
  78 #include "AliRun.h"
  79 #include "AliMUONDigitizerv1.h"
  80
  81
  82 // Defaults parameters for Z positions of chambers
  83 // taken from values for "stations" in AliMUON::AliMUON
  84 //     const Float_t zch[7]={528, 690., 975., 1249., 1449., 1610, 1710.};
  85 // and from array "dstation" in AliMUONv1::CreateGeometry
  86 //          Float_t dstation[5]={20., 20., 20, 20., 20.};
  87 //     for tracking chambers,
  88 //          according to (Z1 = zch - dstation) and  (Z2 = zch + dstation)
  89 //          for the first and second chambers in the station, respectively,
  90 // and from "DTPLANES" in AliMUONv1::CreateGeometry
  91 //           const Float_t DTPLANES = 15.;
  92 //     for trigger chambers,
  93 //          according to (Z1 = zch) and  (Z2 = zch + DTPLANES)
  94 //          for the first and second chambers in the station, respectively
  95
  96 ClassImp(AliMUON)
  97 //__________________________________________________________________
  98 AliMUON::AliMUON()
  99 {
 100 // Default Constructor
 101 //
 102     fNCh             = 0;
 103     fNTrackingCh     = 0;
 104     fIshunt          = 0;
 105     fChambers        = 0;
 106     fGeometryBuilders = 0;
 107     fTriggerCircuits = 0;
 108     fAccMin          = 0.;
 109     fAccMax          = 0.;
 110     fAccCut          = kFALSE;
 111     fMerger          = 0;
 112     fFileName        = 0;
 113     fMUONData        = 0;
 114     fSplitLevel      = 0;
 115 }
 116 //__________________________________________________________________
 117 AliMUON::AliMUON(const char *name, const char *title)
 118   : AliDetector(name,title)
 119 {
 120 //Begin_Html
 121 /*
 122 <img src="gif/alimuon.gif">
 123 */
 124 //End_Html
 125   fMUONData  = 0x0;
 126   fSplitLevel= 0;
 127   fIshunt     =  0;
 128
 129   fNCh             = AliMUONConstants::NCh();
 130   fNTrackingCh     = AliMUONConstants::NTrackingCh();
 131
 132   SetMarkerColor(kRed);//
 133 //
 134 // Creating List of Chambers
 135     Int_t ch;
 136     fChambers = new TObjArray(AliMUONConstants::NCh());
 137     fGeometryBuilders = new TObjArray(AliMUONConstants::NCh());
 138
 139     // Loop over stations
 140     for (Int_t st = 0; st < AliMUONConstants::NCh() / 2; st++) {
 141       // Loop over 2 chambers in the station
 142       for (Int_t stCH = 0; stCH < 2; stCH++) {
 143         //
 144         //
 145         //    Default Parameters for Muon Tracking Stations
 146         ch = 2 * st + stCH;
 147         if (ch < AliMUONConstants::NTrackingCh()) {
 148           fChambers->AddAt(new AliMUONChamber(ch),ch);
 149         } else {
 150           fChambers->AddAt(new AliMUONChamberTrigger(ch),ch);
 151         }
 152         AliMUONChamber* chamber = (AliMUONChamber*) fChambers->At(ch);
 153         //chamber->SetGid(0);
 154         // Default values for Z of chambers
 155         chamber->SetZ(AliMUONConstants::DefaultChamberZ(ch));
 156         //
 157         chamber->InitGeo(AliMUONConstants::DefaultChamberZ(ch));
 158         //          Set chamber inner and outer radius to default
 159         chamber->SetRInner(AliMUONConstants::Dmin(st)/2);
 160         chamber->SetROuter(AliMUONConstants::Dmax(st)/2);
 161         //
 162       } // Chamber stCH (0, 1) in
 163     }     // Station st (0...)
 164
 165     // Negatives values are ignored by geant3 CONS200 in the calculation of the tracking parameters
 166     fMaxStepGas=0.1;
 167     fMaxStepAlu=0.1;
 168     fMaxDestepGas=-1;
 169     fMaxDestepAlu=-1;
 170
 171     fMaxIterPad   = 0;
 172     fCurIterPad   = 0;
 173
 174     fAccMin          = 0.;
 175     fAccMax          = 0.;
 176     fAccCut          = kFALSE;
 177
 178     // cp new design of AliMUONTriggerDecision
 179     fTriggerCircuits = new TObjArray(AliMUONConstants::NTriggerCircuit());
 180     for (Int_t circ=0; circ<AliMUONConstants::NTriggerCircuit(); circ++) {
 181       fTriggerCircuits->AddAt(new AliMUONTriggerCircuit(),circ);
 182     }
 183     fMerger = 0;
 184 }
 185 //____________________________________________________________________
 186 AliMUON::AliMUON(const AliMUON& rMUON):AliDetector(rMUON)
 187 {
 188 // Dummy copy constructor
 189     ;
 190
 191 }
 192 //____________________________________________________________________
 193 AliMUON::~AliMUON()
 194 {
 195 // Destructor
 196   if(fDebug) printf("%s: Calling AliMUON destructor !!!\n",ClassName());
 197   fIshunt  = 0;
 198   if (fMerger) delete fMerger;
 199
 200   if (fGeometryBuilders){
 201     fGeometryBuilders->Delete();
 202     delete fGeometryBuilders;
 203   }
 204   if (fChambers){
 205     fChambers->Delete();
 206     delete fChambers;
 207   }
 208   if (fTriggerCircuits){
 209     fTriggerCircuits->Delete();
 210     delete fTriggerCircuits;
 211   }
 212   delete fMUONData;
 213 }
 214 //_____________________________________________________________________________
 215 void AliMUON::AddGeometryBuilder(AliMUONVGeometryBuilder* geomBuilder)
 216 {
 217 // Adds the geometry builder to the list
 218 // ---
 219
 220   fGeometryBuilders->Add(geomBuilder);
 221 }
 222 //____________________________________________________________________
 223 void AliMUON::BuildGeometry()
 224 {
 225 // Geometry for event display
 226   for (Int_t i=0; i<7; i++) {
 227     for (Int_t j=0; j<2; j++) {
 228       Int_t id=2*i+j+1;
 229       this->Chamber(id-1).SegmentationModel(1)->Draw("eventdisplay");
 230     }
 231   }
 232 }
 233 //___________________________________________________________________
 234 Int_t AliMUON::DistancetoPrimitive(Int_t , Int_t )
 235 {
 236   return 9999;
 237 }
 238 //__________________________________________________________________
 239 void  AliMUON::SetTreeAddress()
 240 {
 241   GetMUONData()->SetLoader(fLoader);
 242   //  GetMUONData()->MakeBranch("D,RC");
 243   //  GetMUONData()->SetTreeAddress("H,D,RC");
 244   GetMUONData()->SetTreeAddress("H");
 245
 246   if (fHits !=  GetMUONData()->Hits())  {
 247     if (gAlice->GetMCApp())
 248       gAlice->GetMCApp()->AddHitList (fHits); // For purifyKine, only necessary when Hit list is created in AliMUONData
 249   }
 250   fHits = GetMUONData()->Hits(); // Added by Ivana to use the methods FisrtHit, NextHit of AliDetector
 251
 252 }
 253
 254 //____________________________________________________________________
 255 void AliMUON::SetPadSize(Int_t id, Int_t isec, Float_t p1, Float_t p2)
 256 {
 257 // Set the pad size for chamber id and cathode isec
 258     Int_t i=2*(id-1);
 259     ((AliMUONChamber*) fChambers->At(i))  ->SetPadSize(isec,p1,p2);
 260     ((AliMUONChamber*) fChambers->At(i+1))->SetPadSize(isec,p1,p2);
 261 }
 262
 263 //___________________________________________
 264 void AliMUON::SetChambersZ(const Float_t *Z)
 265 {
 266   // Set Z values for all chambers (tracking and trigger)
 267   // from the array pointed to by "Z"
 268     for (Int_t ch = 0; ch < AliMUONConstants::NCh(); ch++)
 269         ((AliMUONChamber*) fChambers->At(ch))->SetZ(Z[ch]);
 270     return;
 271 }
 272 //_________________________________________________________________
 273 void AliMUON::SetChambersZToDefault()
 274 {
 275   // Set Z values for all chambers (tracking and trigger)
 276   // to default values
 277   SetChambersZ(AliMUONConstants::DefaultChamberZ());
 278   return;
 279 }
 280 //_________________________________________________________________
 281 void AliMUON::SetChargeSlope(Int_t id, Float_t p1)
 282 {
 283 // Set the inverse charge slope for chamber id
 284     Int_t i=2*(id-1);    //PH    ((AliMUONChamber*) (*fChambers)[i])->SetSigmaIntegration(p1);
 285     //PH    ((AliMUONChamber*) (*fChambers)[i+1])->SetSigmaIntegration(p1);
 286     ((AliMUONChamber*) fChambers->At(i))->SetChargeSlope(p1);
 287     ((AliMUONChamber*) fChambers->At(i+1))->SetChargeSlope(p1);
 288 }
 289 //__________________________________________________________________
 290 void AliMUON::SetChargeSpread(Int_t id, Float_t p1, Float_t p2)
 291 {
 292 // Set sigma of charge spread for chamber id
 293     Int_t i=2*(id-1);
 294     ((AliMUONChamber*) fChambers->At(i))->SetChargeSpread(p1,p2);
 295     ((AliMUONChamber*) fChambers->At(i+1))->SetChargeSpread(p1,p2);
 296 }
 297 //___________________________________________________________________
 298 void AliMUON::SetSigmaIntegration(Int_t id, Float_t p1)
 299 {
 300 // Set integration limits for charge spread
 301     Int_t i=2*(id-1);
 302     ((AliMUONChamber*) fChambers->At(i))->SetSigmaIntegration(p1);
 303     ((AliMUONChamber*) fChambers->At(i+1))->SetSigmaIntegration(p1);
 304 }
 305
 306 //__________________________________________________________________
 307 void AliMUON::SetMaxAdc(Int_t id, Int_t p1)
 308 {
 309 // Set maximum number for ADCcounts (saturation)
 310     Int_t i=2*(id-1);
 311     ((AliMUONChamber*) fChambers->At(i))->SetMaxAdc(p1);
 312     ((AliMUONChamber*) fChambers->At(i+1))->SetMaxAdc(p1);
 313 }
 314 //__________________________________________________________________
 315 void AliMUON::SetMaxStepGas(Float_t p1)
 316 {
 317 // Set stepsize in gas
 318   fMaxStepGas=p1;
 319 }
 320 //__________________________________________________________________
 321 void AliMUON::SetMaxStepAlu(Float_t p1)
 322 {
 323 // Set step size in Alu
 324     fMaxStepAlu=p1;
 325 }
 326 //__________________________________________________________________
 327 void AliMUON::SetMaxDestepGas(Float_t p1)
 328 {
 329 // Set maximum step size in Gas
 330     fMaxDestepGas=p1;
 331 }
 332 //__________________________________________________________________
 333 void AliMUON::SetMaxDestepAlu(Float_t p1)
 334 {
 335 // Set maximum step size in Alu
 336   fMaxDestepAlu=p1;
 337 }
 338 //___________________________________________________________________
 339 void AliMUON::SetAcceptance(Bool_t acc, Float_t angmin, Float_t angmax)
 340 {
 341 // Set acceptance cuts
 342   fAccCut=acc;
 343   fAccMin=angmin*TMath::Pi()/180;
 344   fAccMax=angmax*TMath::Pi()/180;
 345   Int_t ch;
 346   if (acc) {
 347     for (Int_t st = 0; st < AliMUONConstants::NCh() / 2; st++) {
 348       // Loop over 2 chambers in the station
 349       for (Int_t stCH = 0; stCH < 2; stCH++) {
 350         ch = 2 * st + stCH;
 351         //         Set chamber inner and outer radius according to acceptance cuts
 352         Chamber(ch).SetRInner(TMath::Abs(AliMUONConstants::DefaultChamberZ(ch)*TMath::Tan(fAccMin)));
 353         Chamber(ch).SetROuter(TMath::Abs(AliMUONConstants::DefaultChamberZ(ch)*TMath::Tan(fAccMax)));
 354       } // chamber loop
 355     } // station loop
 356   }
 357 }
 358
 359 //____________________________________________________________________
 360 Float_t  AliMUON::GetMaxStepGas() const
 361 {
 362 // Return stepsize in gas
 363
 364   return fMaxStepGas;
 365 }
 366
 367 //____________________________________________________________________
 368 Float_t  AliMUON::GetMaxStepAlu() const
 369 {
 370 // Return step size in Alu
 371
 372   return fMaxStepAlu;
 373 }
 374
 375 //____________________________________________________________________
 376 Float_t  AliMUON::GetMaxDestepGas() const
 377 {
 378 // Return maximum step size in Gas
 379
 380   return fMaxDestepGas;
 381 }
 382
 383 //____________________________________________________________________
 384 Float_t  AliMUON::GetMaxDestepAlu() const
 385 {
 386 // Return maximum step size in Gas
 387
 388   return fMaxDestepAlu;
 389 }
 390
 391 //____________________________________________________________________
 392 void   AliMUON::SetSegmentationModel(Int_t id, Int_t isec, AliSegmentation *segmentation)
 393 {
 394 // Set the segmentation for chamber id cathode isec
 395     ((AliMUONChamber*) fChambers->At(id))->SetSegmentationModel(isec, segmentation);
 396
 397 }
 398 //____________________________________________________________________
 399 void   AliMUON::SetResponseModel(Int_t id, AliMUONResponse *response)
 400 {
 401 // Set the response for chamber id
 402     ((AliMUONChamber*) fChambers->At(id))->SetResponseModel(response);
 403 }
 404 //____________________________________________________________________
 405 void   AliMUON::SetReconstructionModel(Int_t id, AliMUONClusterFinderVS *reconst)
 406 {
 407 // Set ClusterFinder for chamber id
 408     ((AliMUONChamber*) fChambers->At(id))->SetReconstructionModel(reconst);
 409 }
 410 //____________________________________________________________________
 411 void   AliMUON::SetNsec(Int_t id, Int_t nsec)
 412 {
 413 // Set number of segmented cathods for chamber id
 414     ((AliMUONChamber*) fChambers->At(id))->SetNsec(nsec);
 415 }
 416 //____________________________________________________________________
 417 AliDigitizer* AliMUON::CreateDigitizer(AliRunDigitizer* manager) const
 418 {
 419   return new AliMUONDigitizerv1(manager);
 420 }
 421 //_____________________________________________________________________
 422 void AliMUON::SDigits2Digits()
 423 {
 424
 425 // write TreeD here
 426
 427     if (!fMerger) {
 428       if (gAlice->GetDebug()>0) {
 429         cerr<<"AliMUON::SDigits2Digits: create default AliMUONMerger "<<endl;
 430         cerr<<" no merging, just digitization of 1 event will be done"<<endl;
 431       }
 432       fMerger = new AliMUONMerger();
 433     }
 434     fMerger->Init();
 435     fMerger->Digitise();
 436     char hname[30];
 437     //    sprintf(hname,"TreeD%d",fLoader->GetHeader()->GetEvent());
 438     fLoader->TreeD()->Write(hname,TObject::kOverwrite);
 439     fLoader->TreeD()->Reset();
 440 }
 441
 442 //_______________________________________________________________________
 443 AliLoader* AliMUON::MakeLoader(const char* topfoldername)
 444 {
 445 //builds standard getter (AliLoader type)
 446 //if detector wants to use castomized getter, it must overload this method
 447
 448  if (GetDebug())
 449    Info("MakeLoader",
 450         "Creating standard getter for detector %s. Top folder is %s.",
 451          GetName(),topfoldername);
 452  fLoader   = new AliLoader(GetName(),topfoldername);
 453  fMUONData = new AliMUONData(fLoader,GetName(),GetName());
 454  fMUONData->SetSplitLevel(fSplitLevel);
 455  return fLoader;
 456 }
 457
 458 //_______________________________________________________________________
 459 void AliMUON::Trigger(Int_t /*nev*/){
 460 // call the Trigger Algorithm and fill TreeR
 461
 462   Int_t singlePlus[3]  = {0,0,0};
 463   Int_t singleMinus[3] = {0,0,0};
 464   Int_t singleUndef[3] = {0,0,0};
 465   Int_t pairUnlike[3]  = {0,0,0};
 466   Int_t pairLike[3]    = {0,0,0};
 467
 468   ResetTrigger();
 469   AliMUONTriggerDecision* decision= new AliMUONTriggerDecision(fLoader,1);
 470   decision->Trigger();
 471   decision->GetGlobalTrigger(singlePlus, singleMinus, singleUndef,
 472                              pairUnlike, pairLike);
 473
 474   // add a local trigger in the list
 475   GetMUONData()->AddGlobalTrigger(singlePlus, singleMinus, singleUndef, pairUnlike, pairLike);
 476   Int_t i;
 477
 478   for (Int_t icirc=0; icirc<AliMUONConstants::NTriggerCircuit(); icirc++) {
 479     if(decision->GetITrigger(icirc)==1) {
 480       Int_t localtr[7]={0,0,0,0,0,0,0};
 481       Int_t loLpt[2]={0,0}; Int_t loHpt[2]={0,0}; Int_t loApt[2]={0,0};
 482       decision->GetLutOutput(icirc, loLpt, loHpt, loApt);
 483       localtr[0] = icirc;
 484       localtr[1] = decision->GetStripX11(icirc);
 485       localtr[2] = decision->GetDev(icirc);
 486       localtr[3] = decision->GetStripY11(icirc);
 487       for (i=0; i<2; i++) {    // convert the Lut output in 1 digit
 488         localtr[4] = localtr[4]+Int_t(loLpt[i]*TMath::Power(2,i));
 489         localtr[5] = localtr[5]+Int_t(loHpt[i]*TMath::Power(2,i));
 490         localtr[6] = localtr[6]+Int_t(loApt[i]*TMath::Power(2,i));
 491       }
 492       GetMUONData()->AddLocalTrigger(localtr);  // add a local trigger in the list
 493     }
 494   }
 495
 496   delete decision;
 497
 498   //  fLoader->TreeR()->Fill();
 499   GetMUONData()->Fill("GLT"); //Filling Global and Local Trigger GLT
 500   //  char hname[30];
 501   //  sprintf(hname,"TreeR%d",nev);
 502   //  fLoader->TreeR()->Write(hname,TObject::kOverwrite);
 503     //  fLoader->TreeR()->Reset();
 504   fLoader->WriteRecPoints("OVERWRITE");
 505
 506   //  printf("\n End of trigger for event %d\n", nev);
 507 }
 508
 509 //____________________________________________________________________
 510 void AliMUON::Digits2Reco()
 511 {
 512   FindClusters();
 513   Int_t nev = gAlice->GetHeader()->GetEvent();
 514   GetMUONData()->Fill("RC"); //Filling Reconstructed Cluster
 515   fLoader->WriteRecPoints("OVERWRITE");
 516   GetMUONData()->ResetRawClusters();
 517   Info("Digits2Reco","End of cluster finding for event %d", nev);
 518 }
 519 //____________________________________________________________________
 520 void AliMUON::FindClusters()
 521 {
 522 //
 523 //  Perform cluster finding
 524 //
 525     TClonesArray *dig1, *dig2;
 526     Int_t ndig, k;
 527     dig1 = new TClonesArray("AliMUONDigit",1000);
 528     dig2 = new TClonesArray("AliMUONDigit",1000);
 529     AliMUONDigit *digit;
 530 // Loop on chambers and on cathode planes
 531 //
 532     ResetRawClusters();
 533     TClonesArray * muonDigits;
 534
 535     for (Int_t ich = 0; ich < 10; ich++) {
 536       //PH      AliMUONChamber* iChamber = (AliMUONChamber*) (*fChambers)[ich];
 537         AliMUONChamber* iChamber = (AliMUONChamber*) fChambers->At(ich);
 538         AliMUONClusterFinderVS* rec = iChamber->ReconstructionModel();
 539
 540         ResetDigits();
 541         GetMUONData()->GetCathode(0);
 542         //TClonesArray *
 543         muonDigits = GetMUONData()->Digits(ich);
 544         ndig=muonDigits->GetEntriesFast();
 545         if(fDebug)
 546         printf("\n 1 Found %d digits in %p chamber %d", ndig, muonDigits,ich);
 547         TClonesArray &lhits1 = *dig1;
 548         Int_t n = 0;
 549         for (k = 0; k < ndig; k++) {
 550             digit = (AliMUONDigit*) muonDigits->UncheckedAt(k);
 551             if (rec->TestTrack(digit->Track(0)))
 552                 new(lhits1[n++]) AliMUONDigit(*digit);
 553         }
 554         GetMUONData()->ResetDigits();
 555         GetMUONData()->GetCathode(1);
 556         muonDigits =  GetMUONData()->Digits(ich);
 557         ndig=muonDigits->GetEntriesFast();
 558         if(fDebug)
 559         printf("\n 2 Found %d digits in %p %d", ndig, muonDigits, ich);
 560         TClonesArray &lhits2 = *dig2;
 561         n=0;
 562
 563         for (k=0; k<ndig; k++) {
 564             digit= (AliMUONDigit*) muonDigits->UncheckedAt(k);
 565             if (rec->TestTrack(digit->Track(0)))
 566             new(lhits2[n++]) AliMUONDigit(*digit);
 567         }
 568
 569         if (rec) {
 570             AliMUONClusterInput::Instance()->SetDigits(ich, dig1, dig2);
 571             rec->FindRawClusters();
 572         }
 573         dig1->Delete();
 574         dig2->Delete();
 575     } // for ich
 576     delete dig1;
 577     delete dig2;
 578 }
 579 //______________________________________________________________________
 580 #ifdef never
 581 void AliMUON::Streamer(TBuffer &R__b)_
 582 {
 583    // Stream an object of class AliMUON.
 584       AliMUONChamber        *iChamber;
 585       AliMUONTriggerCircuit *iTriggerCircuit;
 586       AliSegmentation       *segmentation;
 587       AliMUONResponse       *response;
 588       TClonesArray          *digitsaddress;
 589       TClonesArray          *rawcladdress;
 590       Int_t i;
 591       if (R__b.IsReading()) {
 592           Version_t R__v = R__b.ReadVersion(); if (R__v) { }
 593           AliDetector::Streamer(R__b);
 594           R__b >> fNPadHits;
 595           R__b >> fPadHits; // diff
 596           R__b >> fNLocalTrigger;
 597           R__b >> fLocalTrigger;
 598           R__b >> fNGlobalTrigger;
 599           R__b >> fGlobalTrigger;
 600           R__b >> fDchambers;
 601           R__b >> fRawClusters;
 602           R__b.ReadArray(fNdch);
 603           R__b.ReadArray(fNrawch);
 604           R__b >> fAccCut;
 605           R__b >> fAccMin;
 606           R__b >> fAccMax;
 607           R__b >> fChambers;
 608           R__b >> fTriggerCircuits;
 609           for (i =0; i<AliMUONConstants::NTriggerCircuit(); i++) {
 610               iTriggerCircuit=(AliMUONTriggerCircuit*) (*fTriggerCircuits)[i];
 611               iTriggerCircuit->Streamer(R__b);
 612           }
 613 // Stream chamber related information
 614           for (i =0; i<AliMUONConstants::NCh(); i++) {
 615               iChamber=(AliMUONChamber*) (*fChambers)[i];
 616               iChamber->Streamer(R__b);
 617               if (iChamber->Nsec()==1) {
 618                   segmentation=iChamber->SegmentationModel(1);
 619                   if (segmentation)
 620                       segmentation->Streamer(R__b);
 621               } else {
 622                   segmentation=iChamber->SegmentationModel(1);
 623                   if (segmentation)
 624                       segmentation->Streamer(R__b);
 625                   if (segmentation)
 626                       segmentation=iChamber->SegmentationModel(2);
 627                   segmentation->Streamer(R__b);
 628               }
 629               response=iChamber->ResponseModel();
 630               if (response)
 631                   response->Streamer(R__b);
 632               digitsaddress=(TClonesArray*) (*fDchambers)[i];
 633               digitsaddress->Streamer(R__b);
 634               if (i < AliMUONConstants::NTrackingCh()) {
 635                   rawcladdress=(TClonesArray*) (*fRawClusters)[i];
 636                   rawcladdress->Streamer(R__b);
 637               }
 638           }
 639
 640       } else {
 641           R__b.WriteVersion(AliMUON::IsA());
 642           AliDetector::Streamer(R__b);
 643           R__b << fNPadHits;
 644           R__b << fPadHits; // diff
 645           R__b << fNLocalTrigger;
 646           R__b << fLocalTrigger;
 647           R__b << fNGlobalTrigger;
 648           R__b << fGlobalTrigger;
 649           R__b << fDchambers;
 650           R__b << fRawClusters;
 651           R__b.WriteArray(fNdch, AliMUONConstants::NCh());
 652           R__b.WriteArray(fNrawch, AliMUONConstants::NTrackingCh());
 653
 654           R__b << fAccCut;
 655           R__b << fAccMin;
 656           R__b << fAccMax;
 657
 658           R__b << fChambers;
 659           R__b << fTriggerCircuits;
 660           for (i =0; i<AliMUONConstants::NTriggerCircuit(); i++) {
 661               iTriggerCircuit=(AliMUONTriggerCircuit*) (*fTriggerCircuits)[i];
 662               iTriggerCircuit->Streamer(R__b);
 663           }
 664           for (i =0; i<AliMUONConstants::NCh(); i++) {
 665               iChamber=(AliMUONChamber*) (*fChambers)[i];
 666               iChamber->Streamer(R__b);
 667               if (iChamber->Nsec()==1) {
 668                   segmentation=iChamber->SegmentationModel(1);
 669                   if (segmentation)
 670                       segmentation->Streamer(R__b);
 671               } else {
 672                   segmentation=iChamber->SegmentationModel(1);
 673                   if (segmentation)
 674                       segmentation->Streamer(R__b);
 675                   segmentation=iChamber->SegmentationModel(2);
 676                   if (segmentation)
 677                       segmentation->Streamer(R__b);
 678               }
 679               response=iChamber->ResponseModel();
 680               if (response)
 681                   response->Streamer(R__b);
 682               digitsaddress=(TClonesArray*) (*fDchambers)[i];
 683               digitsaddress->Streamer(R__b);
 684               if (i < AliMUONConstants::NTrackingCh()) {
 685                   rawcladdress=(TClonesArray*) (*fRawClusters)[i];
 686                   rawcladdress->Streamer(R__b);
 687               }
 688           }
 689       }
 690 }
 691 #endif
 692 //_______________________________________________________________________
 693 AliMUONPadHit* AliMUON::FirstPad(AliMUONHit*  hit, TClonesArray *clusters)
 694 {
 695 // to be removed
 696     // Initialise the pad iterator
 697     // Return the address of the first padhit for hit
 698     TClonesArray *theClusters = clusters;
 699     Int_t nclust = theClusters->GetEntriesFast();
 700     if (nclust && hit->PHlast() > 0) {
 701         AliMUON::fMaxIterPad=hit->PHlast();
 702         AliMUON::fCurIterPad=hit->PHfirst();
 703         return (AliMUONPadHit*) clusters->UncheckedAt(AliMUON::fCurIterPad-1);
 704     } else {
 705         return 0;
 706     }
 707 }
 708 //_______________________________________________________________________
 709 AliMUONPadHit* AliMUON::NextPad(TClonesArray *clusters)
 710 {
 711   // To be removed
 712 // Get next pad (in iterator)
 713 //
 714     AliMUON::fCurIterPad++;
 715     if (AliMUON::fCurIterPad <= AliMUON::fMaxIterPad) {
 716         return (AliMUONPadHit*) clusters->UncheckedAt(AliMUON::fCurIterPad-1);
 717     } else {
 718         return 0;
 719     }
 720 }
 721 //_______________________________________________________________________
 722
 723 AliMUONRawCluster *AliMUON::RawCluster(Int_t ichamber, Int_t icathod, Int_t icluster)
 724 {
 725 //
 726 //  Return rawcluster (icluster) for chamber ichamber and cathode icathod
 727 //  Obsolete ??
 728     TClonesArray *muonRawCluster  = GetMUONData()->RawClusters(ichamber);
 729     ResetRawClusters();
 730     TTree *treeR = fLoader->TreeR();
 731     Int_t nent=(Int_t)treeR->GetEntries();
 732     treeR->GetEvent(nent-2+icathod-1);
 733     //treeR->GetEvent(icathod);
 734     //Int_t nrawcl = (Int_t)muonRawCluster->GetEntriesFast();
 735
 736     AliMUONRawCluster * mRaw = (AliMUONRawCluster*)muonRawCluster->UncheckedAt(icluster);
 737     //printf("RawCluster _ nent nrawcl icluster mRaw %d %d %d%p\n",nent,nrawcl,icluster,mRaw);
 738
 739     return  mRaw;
 740 }
 741 //________________________________________________________________________
 742 void   AliMUON::SetMerger(AliMUONMerger* merger)
 743 {
 744 // Set pointer to merger
 745     fMerger = merger;
 746 }
 747 //________________________________________________________________________
 748 AliMUONMerger*  AliMUON::Merger()
 749 {
 750 // Return pointer to merger
 751     return fMerger;
 752 }
 753 //________________________________________________________________________
 754 AliMUON& AliMUON::operator = (const AliMUON& /*rhs*/)
 755 {
 756 // copy operator
 757 // dummy version
 758     return *this;
 759 }
 760 //________________________________________________________________________
 761 void AliMUON::Reconstruct() const
 762 {
 763
 764 //  AliLoader* loader = GetLoader();
 765
 766   AliRunLoader* runLoader = fLoader->GetRunLoader();
 767   Int_t nEvents = runLoader->GetNumberOfEvents();
 768
 769 // used local container for each method
 770 // passing fLoader as argument, could be avoided ???
 771   AliMUONEventReconstructor* recoEvent = new AliMUONEventReconstructor(fLoader);
 772   AliMUONData* dataEvent = recoEvent->GetMUONData();
 773
 774   AliMUONClusterReconstructor* recoCluster = new AliMUONClusterReconstructor(fLoader);
 775   AliMUONData* dataCluster = recoCluster->GetMUONData();
 776
 777   AliMUONTriggerDecision* trigDec = new AliMUONTriggerDecision(fLoader);
 778   AliMUONData* dataTrig = trigDec->GetMUONData();
 779
 780
 781   for (Int_t i = 0; i < 10; i++) {
 782     AliMUONClusterFinderVS *RecModel = new AliMUONClusterFinderVS();
 783     RecModel->SetGhostChi2Cut(10);
 784     recoCluster->SetReconstructionModel(i,RecModel);
 785   }
 786
 787   fLoader->LoadDigits("READ");
 788   fLoader->LoadRecPoints("RECREATE");
 789   fLoader->LoadTracks("RECREATE");
 790
 791   //   Loop over events
 792   for(Int_t ievent = 0; ievent < nEvents; ievent++) {
 793     printf("Event %d\n",ievent);
 794     runLoader->GetEvent(ievent);
 795
 796     //----------------------- digit2cluster & Digits2Trigger -------------------
 797     if (!fLoader->TreeR()) fLoader->MakeRecPointsContainer();
 798
 799     // tracking branch
 800     dataCluster->MakeBranch("RC");
 801     dataCluster->SetTreeAddress("D,RC");
 802     recoCluster->Digits2Clusters();
 803     dataCluster->Fill("RC");
 804
 805     // trigger branch
 806     dataTrig->MakeBranch("GLT");
 807     dataTrig->SetTreeAddress("D,GLT");
 808     trigDec->Digits2Trigger();
 809     dataTrig->Fill("GLT");
 810
 811     fLoader->WriteRecPoints("OVERWRITE");
 812
 813     //---------------------------- Track & TriggerTrack ---------------------
 814     if (!fLoader->TreeT()) fLoader->MakeTracksContainer();
 815     // tracking branch
 816     dataEvent->MakeBranch("RT"); //track
 817     dataEvent->SetTreeAddress("RT");
 818     recoEvent->EventReconstruct();
 819     dataEvent->Fill("RT");
 820
 821     // trigger branch
 822     dataEvent->MakeBranch("RL"); //trigger track
 823     dataEvent->SetTreeAddress("RL");
 824     recoEvent->EventReconstructTrigger();
 825     dataEvent->Fill("RL");
 826
 827     fLoader->WriteTracks("OVERWRITE");
 828
 829     //--------------------------- Resetting branches -----------------------
 830     dataCluster->ResetDigits();
 831     dataCluster->ResetRawClusters();
 832
 833     dataTrig->ResetDigits();
 834     dataTrig->ResetTrigger();
 835
 836     dataEvent->ResetRawClusters();
 837     dataEvent->ResetTrigger();
 838     dataEvent->ResetRecTracks();
 839     dataEvent->ResetRecTriggerTracks();
 840
 841   }
 842   fLoader->UnloadDigits();
 843   fLoader->UnloadRecPoints();
 844   fLoader->UnloadTracks();
 845
 846   delete recoCluster;
 847   delete recoEvent;
 848   delete trigDec;
 849 }
 850 //________________________________________________________________________
 851 void AliMUON::FillESD(AliESD* event) const
 852 {
 853
 854   TClonesArray* recTracksArray;
 855   TClonesArray* recTrigTracksArray;
 856
 857   //YS AliLoader* loader = GetLoader();
 858   AliRunLoader* runLoader = fLoader->GetRunLoader(); //YS loader->GetRunLoader();
 859   fLoader->LoadTracks("READ"); //YS
 860
 861
 862   // declaration
 863   Int_t ievent;
 864   Int_t ntrackhits;
 865   Double_t fitfmin;
 866
 867
 868   Double_t bendingSlope, nonBendingSlope, fInverseBendingMomentum;
 869   Double_t fXRec, fYRec, fZRec;
 870
 871   Float_t  x11, y11, thetaX,thetaY ;
 872
 873   //YS Int_t nEvents = runLoader->GetNumberOfEvents();
 874
 875   // setting pointer for tracks, triggertracks& trackparam at vertex
 876   AliMUONTrack * rectrack;
 877   AliMUONTriggerTrack * rectriggertrack;
 878   AliMUONTrackParam *trackParam;
 879
 880   ievent = runLoader->GetEventNumber() ; //YS, seems not to be implemented yet (Ch. F)
 881   //YS for (ievent = 0; ievent < nEvents; ievent++) {
 882     runLoader->GetEvent(ievent);
 883
 884     // setting ESD MUON class
 885     AliESDMuonTrack* ESDTrack = new  AliESDMuonTrack() ;
 886
 887     // -------------------- tracks-------------
 888     fMUONData->SetTreeAddress("RT");
 889     fMUONData->GetRecTracks();
 890     recTracksArray = fMUONData->RecTracks();
 891
 892     Int_t nrectracks = (Int_t) recTracksArray->GetEntriesFast(); //
 893
 894     //-------------------- trigger tracks-------------
 895     Bool_t ifirst = true;
 896     Long_t trigPat = 0;
 897     fMUONData->SetTreeAddress("RL");
 898     fMUONData->GetRecTriggerTracks();
 899     recTrigTracksArray = fMUONData->RecTriggerTracks();
 900
 901     Int_t ntrectracks = (Int_t) recTrigTracksArray->GetEntriesFast(); //YS
 902
 903     //printf(">>> Event %d Number of Recconstructed tracks %d \n",ievent, nrectracks);
 904
 905     Int_t ntracks = TMath::Max(nrectracks, ntrectracks);
 906
 907     // loop over tracks
 908     for (Int_t irectracks = 0; irectracks <  ntracks;  irectracks++) {
 909
 910     // -------------------- tracks-------------
 911     // not the best way to do, avoid creating two "branches" in fMuonTracks
 912       if (irectracks < nrectracks) {
 913         rectrack = (AliMUONTrack*) recTracksArray->At(irectracks);
 914
 915         trackParam = rectrack->GetTrackParamAtVertex();
 916
 917         bendingSlope            = trackParam->GetBendingSlope();
 918         nonBendingSlope         = trackParam->GetNonBendingSlope();
 919         fInverseBendingMomentum = trackParam->GetInverseBendingMomentum();
 920         fXRec  = trackParam->GetNonBendingCoor();
 921         fYRec  = trackParam->GetBendingCoor();
 922         fZRec  = trackParam->GetZ();
 923
 924         ntrackhits = rectrack->GetNTrackHits();
 925         fitfmin = rectrack->GetFitFMin();
 926
 927         // setting data member of ESD MUON
 928         ESDTrack->SetInverseBendingMomentum(fInverseBendingMomentum);
 929         ESDTrack->SetThetaX(TMath::ATan(nonBendingSlope));
 930         ESDTrack->SetThetaY(TMath::ATan(bendingSlope));
 931         ESDTrack->SetZ(fZRec);
 932         ESDTrack->SetBendingCoor(fYRec);
 933         ESDTrack->SetNonBendingCoor(fXRec);
 934         ESDTrack->SetChi2(fitfmin);
 935         ESDTrack->SetNHit(ntrackhits);
 936       }
 937     // -------------------- trigger tracks-------------
 938       if (irectracks < ntrectracks) {
 939         rectriggertrack = (AliMUONTriggerTrack*) recTrigTracksArray->At(irectracks);
 940         if (ifirst) {
 941           ifirst = false;
 942           trigPat = rectriggertrack->GetGTPattern();
 943         }
 944
 945         x11 = rectriggertrack->GetY11();
 946         y11 = rectriggertrack->GetY11();
 947         thetaX = rectriggertrack->GetThetax();
 948         thetaY = rectriggertrack->GetThetay();
 949
 950         // setting data member of ESD MUON trigger
 951         ESDTrack->SetThetaX11(thetaX);
 952         ESDTrack->SetThetaY11(thetaY);
 953         ESDTrack->SetX11(x11);
 954         ESDTrack->SetY11(y11);
 955       }
 956       // storing ESD MUON Track into ESD Event & reset muondata
 957       if (ntrectracks+ntrectracks != 0)
 958         event->AddMuonTrack(ESDTrack);
 959     } // end loop tracks
 960
 961     //    if (ntrectracks+ntrectracks != 0)
 962       //  event->SetTrigger(trigPat);
 963       // Waiting for STEER update
 964     fMUONData->ResetRecTracks();
 965     fMUONData->ResetRecTriggerTracks();
 966
 967     //} // end loop on event
 968     fLoader->UnloadTracks();
 969 }
 970