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