New Stepmanager function.
authormartinez <martinez@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 9 Oct 2003 13:43:57 +0000 (13:43 +0000)
committermartinez <martinez@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 9 Oct 2003 13:43:57 +0000 (13:43 +0000)
(Khalil Boudjemline Gines Martinez)

New coordinate system in Slat segmentation
(Eric Dumonteil)

Config file option param operational

13 files changed:
MUON/AliMUON.cxx
MUON/AliMUONData.cxx
MUON/AliMUONData.h
MUON/AliMUONDigitizerv1.cxx
MUON/AliMUONDisplay.cxx
MUON/AliMUONEventReconstructor.cxx
MUON/AliMUONHit.cxx
MUON/AliMUONHit.h
MUON/AliMUONSegmentationSlat.cxx
MUON/AliMUONv1.cxx
MUON/AliMUONv1.h
MUON/Config_MUON_test.C
MUON/MUONrawclusters.C

index d34c95898abb1c8a13ae2a284ae81e98492d86b7..addbdb55ee8ed01db6e05029066819454cf81f21 100644 (file)
@@ -147,8 +147,9 @@ AliMUON::AliMUON(const char *name, const char *title)
       } // Chamber stCH (0, 1) in 
     }     // Station st (0...)
     
-    fMaxStepGas=0.01; 
-    fMaxStepAlu=0.1; 
+    // Negatives values are ignored by geant3 CONS200 in the calculation of the tracking parameters
+    fMaxStepGas=0.1; 
+    fMaxStepAlu=0.1;  
     fMaxDestepGas=-1;
     fMaxDestepAlu=-1;
     
index e3b70b483a43dedeb1b29b2216431e5c3211f5f5..fa6a74ebd3a2e119b297540d142607cee475edf0 100644 (file)
@@ -150,6 +150,20 @@ void AliMUONData::AddHit(Int_t fIshunt, Int_t track, Int_t iChamber,
                                  phi, length, destep);
 }
 //____________________________________________________________________________
+void AliMUONData::AddHit(Int_t fIshunt, Int_t track, Int_t iChamber, 
+                        Int_t idpart, Float_t X, Float_t Y, Float_t Z, 
+                        Float_t tof, Float_t momentum, Float_t theta, 
+                        Float_t phi, Float_t length, Float_t destep,
+                        Float_t Xref,Float_t Yref,Float_t Zref)
+{
+  TClonesArray &lhits = *fHits;
+  new(lhits[fNhits++]) AliMUONHit(fIshunt, track, iChamber, 
+                                 idpart, X, Y, Z, 
+                                 tof, momentum, theta, 
+                                 phi, length, destep,
+                                 Xref,Yref,Zref);
+}
+//____________________________________________________________________________
 void AliMUONData::AddLocalTrigger(Int_t *localtr)
 {
   // add a MUON Local Trigger to the list
index f54c0cdf00f7721cb7c695797a3f67f5b113f9fd..abbb0a94fc7f38d73d7be728f69537ca70aac9ae 100644 (file)
@@ -1,7 +1,7 @@
 #ifndef ALIMUONDATA_H
 #define ALIMUONDATA_H
 //
-/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+/*  Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
  * See cxx source for full Copyright notice                               */
 
 /* $Id$ */
@@ -42,6 +42,11 @@ class AliMUONData : public TNamed {
                          Int_t idpart, Float_t X, Float_t Y, Float_t Z, 
                          Float_t tof, Float_t momentum, Float_t theta, 
                          Float_t phi, Float_t length, Float_t destep);
+    virtual void   AddHit(Int_t fIshunt, Int_t track, Int_t iChamber, 
+                         Int_t idpart, Float_t X, Float_t Y, Float_t Z, 
+                         Float_t tof, Float_t momentum, Float_t theta, 
+                         Float_t phi, Float_t length, Float_t destep, 
+                         Float_t Xref,Float_t Yref,Float_t Zref);
     virtual void   AddGlobalTrigger(Int_t *singlePlus, Int_t *singleMinus,
                                    Int_t *singleUndef, Int_t *pairUnlike, 
                                    Int_t *pairLike);
index ceafc77e6a78f12d2ede99546330209fdbbfb36b..882c55b06ac9e97230a6de8b2ed7aff4d287730f 100644 (file)
@@ -77,7 +77,25 @@ Bool_t AliMUONDigitizerv1::ExistTransientDigit(AliMUONTransientDigit * mTD)
 //------------------------------------------------------------------------
 Bool_t AliMUONDigitizerv1::Init()
 {
+
 // Initialization 
+  if (GetDebug()>2) Info("Init","AliMUONDigitizerv1::Init() starts");
+
+  //Loaders (We assume input0 to be the output too)
+  AliRunLoader * runloader;    // Input   loader
+  AliLoader    * gime; 
+  // Getting runloader
+  runloader    = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
+  if (runloader == 0x0) {
+    Error("Init","RunLoader is not in input file 0");
+    return kFALSE;  // RunDigitizer is not working.  
+  }
+  // Getting MUONloader
+  gime        = runloader->GetLoader("MUONLoader");
+  gime->LoadHits("READ");
+  gime->LoadDigits("RECREATE");
   return kTRUE;
 }
 
index 5c8a018d7e4dde2ce8833c363009a897d882bfd5..1dc4250c045bbb1bd9c05478824272f886940586 100644 (file)
@@ -622,6 +622,7 @@ void AliMUONDisplay::DrawTitle(Option_t *option)
     Float_t dy   = ymax-ymin;
     
 
+    AliRunLoader * RunLoader = AliRunLoader::GetRunLoader("MUONFolder");
 
     if (strlen(option) == 0) {
        TPaveText *title = new TPaveText(xmin +0.01*dx, ymax-0.09*dy, xmin +0.5*dx, ymax-0.01*dy);
@@ -632,7 +633,7 @@ void AliMUONDisplay::DrawTitle(Option_t *option)
        title->Draw();
        char ptitle[100];
        sprintf(ptitle, "Alice event:%d Run:%d Chamber:%d Cathode:%d",
-               gAlice->GetHeader()->GetEvent(),
+               RunLoader->GetEventNumber(),
                gAlice->GetHeader()->GetRun(),
                fChamber,
                fCathode);
@@ -1119,18 +1120,20 @@ void AliMUONDisplay::SetView(Float_t theta, Float_t phi, Float_t psi)
 //_____________________________________________________________________________
 void AliMUONDisplay::ShowNextEvent(Int_t delta)
 {
+
+  AliRunLoader * RunLoader = AliRunLoader::GetRunLoader("MUONFolder");
+    
 //  Display (current event_number + delta)
 //    delta =  1  shown next event
 //    delta = -1 show previous event
     if (delta) {
-       gAlice->Clear();
-       Int_t currentEvent = gAlice->GetHeader()->GetEvent();
-       Int_t newEvent     = currentEvent + delta;
-       gAlice->GetEvent(newEvent);
-       fEvent=newEvent;
-       if (!gAlice->TreeD()) return; 
+      //RunLoader->CleanDetectors();
+      //RunLoader->CleanKinematics();
+      Int_t currentEvent = RunLoader->GetEventNumber();
+      Int_t newEvent     = currentEvent + delta;
+      RunLoader->GetEvent(newEvent);
+      fEvent=newEvent;
     }
-    
     LoadDigits(fChamber, fCathode);
     fPad->cd(); 
     Draw();
index af219e38f23c98546b3ecdce6e7b3a7a64ed91e3..c79f995e21d0fd05b6cffd4a1f360f9010b30dda 100644 (file)
@@ -509,13 +509,15 @@ void AliMUONEventReconstructor::AddHitsForRecFromGEANT(TTree *TH)
     Int_t ihit, nhits=0;
       nhits = (Int_t) muondata->Hits()->GetEntriesFast();
       AliMUONHit* mHit=0x0;
+
       for(ihit=0; ihit<nhits; ihit++) {
        mHit = static_cast<AliMUONHit*>(muondata->Hits()->At(ihit));
        Int_t ipart = TMath::Abs ((Int_t) mHit->Particle()); //AZ
        if (NewHitForRecFromGEANT(mHit,track, hit, 1) && ipart == 13
            //if (NewHitForRecFromGEANT(mHit,itrack-1, hit, 1) && ipart == 13 
-           && itrack <= 2) chamBits |= BIT(mHit->Chamber()-1); //AZ - set bit
+           && itrack <= 2 && !BIT(mHit->Chamber()-1)  ) chamBits |= BIT(mHit->Chamber()-1); //AZ - set bit
       }
+
     if (chamBits&3 && chamBits>>2&3 && chamBits>>4&3 && chamBits>>6&3 && 
         chamBits>>8&3 && ((chamBits>>6&3)==3 || (chamBits>>8&3)==3)) 
       fMuons += 1; //AZ
index d59ea595806abbc87350fc0424d3f58802bc9268..00d6d90fde398f3b39756a90bbe38c9df48aa311 100644 (file)
@@ -40,9 +40,14 @@ AliMUONHit::AliMUONHit(Int_t shunt, Int_t track, Int_t *vol, Float_t *hits):
     fPy        = hits[12];
     fPz        = hits[13];
     fAge       = hits[14];
+    fXref      = 0.;
+    fYref      = 0.;
+    fZref      = 0.;
 }
 //___________________________________________
-AliMUONHit::AliMUONHit(Int_t shunt, Int_t track, Int_t iChamber, Int_t idpart, Float_t X, Float_t Y, Float_t Z, Float_t tof, Float_t momentum, Float_t theta, Float_t phi, Float_t length, Float_t destep):
+AliMUONHit::AliMUONHit(Int_t shunt, Int_t track, Int_t iChamber, Int_t idpart, 
+                      Float_t X, Float_t Y, Float_t Z, Float_t tof, Float_t momentum, 
+                      Float_t theta, Float_t phi, Float_t length, Float_t destep):
        AliHit(shunt, track)
 {
 // Constructor
@@ -59,12 +64,34 @@ AliMUONHit::AliMUONHit(Int_t shunt, Int_t track, Int_t iChamber, Int_t idpart, F
     fPHlast    = 0;
     fPTot      = momentum;
     fAge       = tof;
+    fXref      = 0.;
+    fYref      = 0.;
+    fZref      = 0.;
 }
-
-
-
-
-
-
-
+//-----------------------------------------------------------------------------------------------
+AliMUONHit::AliMUONHit(Int_t shunt, Int_t track, Int_t iChamber, Int_t idpart, 
+                      Float_t X, Float_t Y, Float_t Z, Float_t tof, Float_t momentum, 
+                      Float_t theta, Float_t phi, Float_t length, Float_t destep,
+                      Float_t Xref,Float_t Yref,Float_t Zref):
+       AliHit(shunt, track)
+{
+// Constructor
+    fChamber   = iChamber;
+    fParticle  = idpart;
+    fX         = X;
+    fY         = Y;
+    fZ         = Z;
+    fTheta     = theta;
+    fPhi       = phi;
+    fTlength   = length;
+    fEloss     = destep;
+    fPHfirst   = 0;
+    fPHlast    = 0;
+    fPTot      = momentum;
+    fAge       = tof;
+    fXref      = Xref;
+    fYref      = Yref;
+    fZref      = Zref;
+}
+//-----------------------------------------------------------------------------------------------
 
index edd1e4d70db5fe8d61dd3f6c2325b1e22bd0aa1b..44c2cf6bdf0d134d373934fcb54ae4c52d73f6f5 100644 (file)
@@ -14,6 +14,10 @@ class AliMUONHit : public AliHit {
     AliMUONHit() {}
     AliMUONHit(Int_t fIshunt, Int_t track, Int_t *vol, Float_t *hits);
     AliMUONHit(Int_t fIshunt, Int_t track, Int_t iChamber, Int_t idpart, Float_t X, Float_t Y, Float_t Z, Float_t tof, Float_t momentum, Float_t theta, Float_t phi, Float_t length, Float_t destep);
+    AliMUONHit(Int_t fIshunt, Int_t track, Int_t iChamber, Int_t idpart, 
+               Float_t X, Float_t Y, Float_t Z, Float_t tof, Float_t momentum, 
+               Float_t theta, Float_t phi, Float_t length, Float_t destep,
+               Float_t Xref, Float_t Yref, Float_t Zref);
     virtual ~AliMUONHit() {}
     Int_t   Chamber()  {return fChamber;}
     Float_t Particle() {return fParticle;}    
@@ -32,6 +36,11 @@ class AliMUONHit : public AliHit {
     Float_t Cy()       {return fPy/fPTot;}
     Float_t Cz()       {return fPz/fPTot;}
 
+    Float_t Xref()     {return fXref;}
+    Float_t Yref()     {return fYref;}
+    Float_t Zref()     {return fZref;}
+
+
  private:
     Int_t     fChamber;       // Chamber number
     Float_t   fParticle;      // Geant3 particle type
@@ -43,11 +52,16 @@ class AliMUONHit : public AliHit {
     Int_t     fPHfirst;       // first padhit
     Int_t     fPHlast;        // last padhit
 
-    Float_t   fPTot;          // hit local momentum P
+    Float_t   fPTot;          // Local momentum P of the track when entering in the chamber
     Float_t   fPx;            // Px
     Float_t   fPy;            // Py
     Float_t   fPz;            // Pz
     
+    Float_t   fXref;          // X position of hit in the center of the chamber (without angle effect)
+    Float_t   fYref;          // Y position of hit in the center of the chamber (without angle effect)
+    Float_t   fZref;          // Z position of hit in the center of the chamber (without angle effect)
+
+    
     ClassDef(AliMUONHit,1)    //Hit object for MUON
 };
 #endif
index 9add3aeb9ac7df61a93af7ca8c6a22a39b37c93f..9e81968159293b09208003d15428629d5909a52e 100644 (file)
@@ -137,15 +137,17 @@ void AliMUONSegmentationSlat::GlobalToLocal(
 //                                                 positive side is shifted up
 // by half the overlap
     zlocal = z-fChamber->Z();
-    zlocal = (x>0) ? zlocal-2.*fDz : zlocal+2.*fDz;
+
+//     zlocal = (x>0) ? zlocal-2.*fDz : zlocal+2.*fDz;
+    zlocal = (x>0) ? zlocal+2.*fDz : zlocal-2.*fDz;      //Change?
+
 //  Set the signs for the symmetry transformation and transform to first quadrant
     SetSymmetry(x);
     Float_t xabs=TMath::Abs(x);
 
-    Int_t ifirst = (zlocal < Float_t(0))? 0:1;
-//
+
 // Find slat number                      
-    for (i=ifirst; i<fNSlats; i+=2) {
+    for (i=0; i<fNSlats; i+=1) {       //Loop on all slats (longuer but more secure)
        index=i;
        if ((y >= fYPosition[i]-eps) && (y <= fYPosition[i]+fSlatY+eps)) break;
     }
@@ -203,7 +205,7 @@ LocalToGlobal(Int_t islat, Float_t  xlocal, Float_t  ylocal, Float_t  &x, Float_
     x = (xlocal+fXPosition[islat])*fSym;
     y=(ylocal+fYPosition[islat]);
 
-    z = (TMath::Even(islat)) ?     -fDz : fDz ; 
+    z = (TMath::Even(islat)) ?     fDz : -fDz ; //Change for new referential
     z = (x>0)                ? z+2.*fDz : z-2.*fDz ; 
 
     z+=fChamber->Z();
@@ -277,7 +279,7 @@ GetPadC(Int_t ix, Int_t iy, Float_t &x, Float_t &y, Float_t &z)
     x=x*TMath::Sign(1,ix);
 
 // z-position
-    z = (TMath::Even(islat)) ?      -fDz : fDz ; 
+    z = (TMath::Even(islat)) ?      fDz : -fDz ; //Change for new referential
     z = (x>0)                ?  z+2.*fDz : z-2.*fDz ; 
     z += fChamber->Z();
 }
index 2cd88bec14fecdc699ef12623f2b918c117f4414..ea11291fe593f8604923b46332da6f6c1df2c469 100644 (file)
@@ -27,6 +27,7 @@
 #include <TRandom.h> 
 #include <TTUBE.h>
 #include <TVirtualMC.h>
+#include <TParticle.h>
 
 #include "AliCallf77.h"
 #include "AliConst.h" 
@@ -44,21 +45,22 @@ ClassImp(AliMUONv1)
  
 //___________________________________________
 AliMUONv1::AliMUONv1() : AliMUON()
+  ,fTrackMomentum(), fTrackPosition()
 {
 // Constructor
-    fChambers = 0;
-    fStations = 0;
-    fStepManagerVersionOld = kFALSE;
-    fStepManagerVersionNew = kFALSE;
-    fStepManagerVersionTest = kFALSE;
-
-    fStepMaxInActiveGas = 2.0;
-}
-
+    fChambers   = 0;
+    fStations   = 0;
+    fStepManagerVersionOld  = kFALSE;
+    fStepMaxInActiveGas     = 0.6;
+    fStepSum    =  0x0;
+    fDestepSum  =  0x0;
+    fElossRatio =  0x0;
+    fAngleEffect10   = 0x0;
+    fAngleEffectNorma= 0x0;
+} 
 //___________________________________________
 AliMUONv1::AliMUONv1(const char *name, const char *title)
-       : AliMUON(name,title)
+  : AliMUON(name,title), fTrackMomentum(), fTrackPosition()
 {
 // Constructor
     // By default include all stations
@@ -69,10 +71,35 @@ AliMUONv1::AliMUONv1(const char *name, const char *title)
     factory.Build(this, title);
 
     fStepManagerVersionOld = kFALSE;
-    fStepManagerVersionNew = kFALSE;
-    fStepManagerVersionTest = kFALSE;
 
-    fStepMaxInActiveGas = 2.0;
+    fStepMaxInActiveGas = 0.6;
+
+    fStepSum   = new Float_t [AliMUONConstants::NCh()];
+    fDestepSum = new Float_t [AliMUONConstants::NCh()];
+    for (Int_t i=0; i<AliMUONConstants::NCh(); i++) {
+      fStepSum[i] =0.0;
+      fDestepSum[i]=0.0;
+    }
+    // Ratio of particle mean eloss with respect MIP's Khalil Boudjemline, sep 2003, PhD.Thesis and Particle Data Book
+    fElossRatio = new TF1("ElossRatio","[0]+[1]*x+[2]*x*x+[3]*x*x*x+[4]*x*x*x*x",0.5,5.); 
+    fElossRatio->SetParameter(0,1.02138);
+    fElossRatio->SetParameter(1,-9.54149e-02);
+    fElossRatio->SetParameter(2,+7.83433e-02); 
+    fElossRatio->SetParameter(3,-9.98208e-03);
+    fElossRatio->SetParameter(4,+3.83279e-04);
+
+    // Angle effect in tracking chambers at theta =10 degres as a function of ElossRatio (Khalil BOUDJEMLINE sep 2003 Ph.D Thesis) (in micrometers)
+    fAngleEffect10 = new TF1("AngleEffect10","[0]+[1]*x+[2]*x*x",0.5,3.0);
+    fAngleEffect10->SetParameter(0, 1.90691e+02);
+    fAngleEffect10->SetParameter(1,-6.62258e+01);
+    fAngleEffect10->SetParameter(2,+1.28247e+01);
+    // Angle effect: Normalisation form theta=10 degres to theta between 0 and 10 (Khalil BOUDJEMLINE sep 2003 Ph.D Thesis)  
+    // Angle with respect to the wires assuming that chambers are perpendicular to the z axis.
+    fAngleEffectNorma = new TF1("AngleEffectNorma","[0]+[1]*x+[2]*x*x+[3]*x*x*x",0.0,10.0);
+    fAngleEffectNorma->SetParameter(0,4.148);
+    fAngleEffectNorma->SetParameter(1,-6.809e-01);
+    fAngleEffectNorma->SetParameter(2,5.151e-02);
+    fAngleEffectNorma->SetParameter(3,-1.490e-03);
 }
 
 //___________________________________________
@@ -1596,190 +1623,139 @@ void AliMUONv1::Init()
    //cp
 
 }
-//___________________________________________
+
+//_______________________________________________________________________________
+Int_t  AliMUONv1::GetChamberId(Int_t volId) const
+{
+// Check if the volume with specified  volId is a sensitive volume (gas) 
+// of some chamber and returns the chamber number;
+// if not sensitive volume - return 0.
+// ---
+
+  for (Int_t i = 1; i <= AliMUONConstants::NCh(); i++)
+    if (volId==((AliMUONChamber*)(*fChambers)[i-1])->GetGid()) return i;
+
+  return 0;
+}
+//_______________________________________________________________________________
 void AliMUONv1::StepManager()
 {
 if (fStepManagerVersionOld) {
+ if (fStepManagerVersionOld) {
     StepManagerOld();
     return;
   }
- if (fStepManagerVersionNew) {
-    StepManagerNew();
-    return;
-  }
- if (fStepManagerVersionTest) {
-    StepManagerTest();
-    return;
-  }
-
-
-  // Volume id
-  Int_t   copy, id;
-  Int_t   idvol;
-  Int_t   iChamber=0;
-  // Particule id, pos and mom vectors, 
-  // theta, phi angles with respect the normal of the chamber, 
-  // spatial step, delta_energy and time of flight
-  Int_t          ipart;
-  TLorentzVector pos, mom;
-  Float_t        theta, phi, tof;
-  Float_t        destep, step;
-  const  Float_t kBig = 1.e10;
 
   // Only charged tracks
   if( !(gMC->TrackCharge()) ) return; 
-
+  // Only charged tracks
+  
   // Only gas gap inside chamber
   // Tag chambers and record hits when track enters 
-  idvol=-1;
+  Int_t   idvol=-1;
+  Int_t   iChamber=0;
+  Int_t   id=0;
+  Int_t   copy;
+  const  Float_t kBig = 1.e10;
+
   id=gMC->CurrentVolID(copy);
+  // printf("id == %d \n",id);
   for (Int_t i = 1; i <= AliMUONConstants::NCh(); i++) {
     if(id==((AliMUONChamber*)(*fChambers)[i-1])->GetGid()) {
       iChamber = i;
       idvol  = i-1;
     }
   }
-  if (idvol == -1) return;
-
-  //   printf(">>>> This Chamber %d\n",iChamber);
-
-  // record hits when track enters ...
-  if( gMC->IsTrackEntering()) gMC->SetMaxStep(fStepMaxInActiveGas);  
-
-  if (gMC->TrackStep() > 0.) {
-    // Get current particle id (ipart), track position (pos)  and momentum (mom)
-    gMC->TrackPosition(pos);
-    gMC->TrackMomentum(mom);
-    ipart    = gMC->TrackPid();
-    theta    = mom.Theta()*kRaddeg;     // theta of track
-    phi      = mom.Phi()  *kRaddeg;     // phi of the track
-    tof      = gMC->TrackTime();        // Time of flight
-    //
-    // momentum loss and steplength in last step
-    destep = gMC->Edep();
-    step   = gMC->TrackStep();
-
-    //new hit
-    GetMUONData()->AddHit(fIshunt, gAlice->GetCurrentTrackNumber(), iChamber, ipart, 
-                         pos.X(), pos.Y(), pos.Z(), tof, mom.P(), 
-                         theta, phi, step, destep);
-  }
-  // Track left chamber ...
-  if( gMC->IsTrackExiting() || gMC->IsTrackStop() || gMC->IsTrackDisappeared()){
-    gMC->SetMaxStep(kBig);
+  if (idvol == -1) {
+    return;
   }
-}
-
-
-Int_t  AliMUONv1::GetChamberId(Int_t volId) const
-{
-// Check if the volume with specified  volId is a sensitive volume (gas) 
-// of some chamber and returns the chamber number;
-// if not sensitive volume - return 0.
-// ---
 
-  for (Int_t i = 1; i <= AliMUONConstants::NCh(); i++)
-    if (volId==((AliMUONChamber*)(*fChambers)[i-1])->GetGid()) return i;
-
-  return 0;
-}
-//__
-
-
-
-void AliMUONv1::StepManagerTest()
-{
-  return;
-}
-//________________________________________
-void AliMUONv1::StepManagerNew()
-{
-
-//   // Volume id
-//   Int_t   copy, id;
-//   Int_t   idvol;
-//   Int_t   iChamber=0;
-//   // Particule id, pos and mom vectors, 
-//   // theta, phi angles with respect the normal of the chamber, 
-//   // spatial step, delta_energy and time of flight
-//   Int_t          ipart;
-//   TLorentzVector pos, mom;
-//   Float_t        theta, phi, tof;
-//   Float_t        destep, step;
-//   const  Float_t kBig = 1.e10;
-
-//   // Only charged tracks
-//   if( !(gMC->TrackCharge()) ) return; 
-
-//   // Only gas gap inside chamber
-//   // Tag chambers and record hits when track enters 
-//   idvol=-1;
-//   id=gMC->CurrentVolID(copy);
-//   for (Int_t i = 1; i <= AliMUONConstants::NCh(); i++) {
-//     if(id==((AliMUONChamber*)(*fChambers)[i-1])->GetGid()) {
-//       iChamber = i;
-//       idvol  = i-1;
-//     }
-//   }
-//   static Float_t Sstep[20]; // Sum of steps per chamber
-//   // static Float_t Sdestep[20]; // Sum of eloss per chamber
-//   Float_t GAP;
-//   Float_t TEST; 
-
-//   if (idvol == -1) return;
-
-//   //  printf(">>>> This Chamber %d\n",iChamber);
-
-//   // record hits when track enters ...
-//   //if( gMC->IsTrackEntering()) gMC->SetMaxStep(fStepMaxInActiveGas);  
-
-//   if (gMC->TrackStep() > 0.) {
-//     // Get current particle id (ipart), track position (pos)  and momentum (mom)
-//     gMC->TrackPosition(pos);
-//     gMC->TrackMomentum(mom);
-//     ipart    = gMC->TrackPid();         // Particle
-//     theta    = mom.Theta()*kRaddeg;     // theta of track
-//     phi      = mom.Phi()  *kRaddeg;     // phi of the track
-//     tof      = gMC->TrackTime();        // Time of flight
-//     //
-//     // momentum loss and steplength in last step
-//     destep = gMC->Edep();
-//     step   = gMC->TrackStep();
-
-//     Sstep[iChamber]+=step;
-//     // Sdestep[iChamber]+=destep;
+   if( gMC->IsTrackEntering() ) {
+     Float_t theta = fTrackMomentum.Theta();
+     Float_t phi   = fTrackMomentum.Phi();
+     Float_t theta_wires      =  TMath::Abs( TMath::ASin( TMath::Sin(theta) * TMath::Sin(phi) ) );
+     if ( (theta_wires*kRaddeg>=10) ) gMC->SetMaxStep(fStepMaxInActiveGas);
+   }
 
-//   }
+//  if (GetDebug()) {
+//     Float_t z = ( (AliMUONChamber*)(*fChambers)[idvol])->Z() ;
+//      Info("StepManager Step","Active volume found %d chamber %d Z chamber is %f ",idvol,iChamber, z);
+//   }  
+  // Particule id and mass, 
+  Int_t     ipart = gMC->TrackPid();
+  Float_t   mass  = gMC->TrackMass();
+
+  fDestepSum[idvol]+=gMC->Edep();
+  // Get current particle id (ipart), track position (pos)  and momentum (mom)
+  if ( fStepSum[idvol]==0.0 )  gMC->TrackMomentum(fTrackMomentum);
+  fStepSum[idvol]+=gMC->TrackStep();
   
-//   step   =  Sstep[iChamber]; // Total step >= gap
-//   //  destep   =  Sdestep[iChamber]; // Total eloss
-
-
-//   // Track left chamber ...
-//   if( gMC->IsTrackExiting() || gMC->IsTrackStop() || gMC->IsTrackDisappeared()){
-//     gMC->SetMaxStep(kBig);
-
-//     Sstep[iChamber]=0; // Reset for the next event
-//     //Sdestep[iChamber]=0; // Reset for the next event
+//   if (GetDebug()) {
+//     Info("StepManager Step","iChamber %d, Particle %d, theta %f phi %f mass %f StepSum %f eloss %g",
+//      iChamber,ipart, fTrackMomentum.Theta()*kRaddeg, fTrackMomentum.Phi()*kRaddeg, mass, fStepSum[idvol], gMC->Edep());
+//     Info("StepManager Step","Track Momentum %f %f %f", fTrackMomentum.X(), fTrackMomentum.Y(), fTrackMomentum.Z()) ;
+//     gMC->TrackPosition(fTrackPosition);
+//     Info("StepManager Step","Track Position %f %f %f",fTrackPosition.X(),fTrackPosition.Y(),fTrackPosition.Z()) ;
+//   }
 
-//     if (iChamber>=1 && iChamber<=2) GAP=0.4;
-//     if (iChamber>=11 && iChamber<=14) GAP=0.2;
-//     if (iChamber>=3 && iChamber<=10) GAP=0.5;
-   
-//     TF1 *ELOSS1 = new TF1("Gauss1","exp(-((x-4.13727e+01)**2)/(2*1.42223e+01**2))",0,75);
-//     TF1 *ELOSS2 = new TF1("Gauss2","exp(-((x+6.83795e+02)**2)/(2*4.48415e+02**2))",75,350);
-//     TEST=gRandom->Rndm();
-//     if (TEST <=0.89) destep=ELOSS1->GetRandom();
-//     else destep=ELOSS2->GetRandom();
-//     destep*=pow(10,-6)*0.0274;
-//     destep*=GAP/0.5;
-
-//     // One hit per chamber
-//     GetMUONData()->AddHit(fIshunt, gAlice->GetCurrentTrackNumber(), iChamber, ipart, 
-//        pos.X()-(step/2*sin(theta*kDegrad)*cos(phi*kDegrad)), pos.Y()-(step/2*sin(theta*kDegrad)*sin(phi*kDegrad)), pos.Z()-GAP/2, tof, mom.P(),theta, phi, step, destep);
+  // Track left chamber or StepSum larger than fStepMaxInActiveGas
+  if ( gMC->IsTrackExiting() || 
+       gMC->IsTrackStop() || 
+       gMC->IsTrackDisappeared()||
+       (fStepSum[idvol]>fStepMaxInActiveGas) ) {
+    
+    if   ( gMC->IsTrackExiting() || 
+          gMC->IsTrackStop() || 
+          gMC->IsTrackDisappeared() ) gMC->SetMaxStep(kBig);
 
-//   }
+    gMC->TrackPosition(fTrackPosition);
+    Float_t theta = fTrackMomentum.Theta();
+    Float_t phi   = fTrackMomentum.Phi();
+    
+    TLorentzVector BackToWire( fStepSum[idvol]/2.*sin(theta)*cos(phi),
+                              fStepSum[idvol]/2.*sin(theta)*sin(phi),
+                              fStepSum[idvol]/2.*cos(theta),0.0       );
+    //     if (GetDebug()) 
+    //       Info("StepManager Exit","Track Position %f %f %f",fTrackPosition.X(),fTrackPosition.Y(),fTrackPosition.Z()) ;
+    //     if (GetDebug()) 
+    //        Info("StepManager Exit ","Track BackToWire %f %f %f",BackToWire.X(),BackToWire.Y(),BackToWire.Z()) ;
+    fTrackPosition-=BackToWire;
+    
+    //-------------- Angle effect 
+    // Ratio between energy loss of particle and Mip as a function of BetaGamma of particle (Energy/Mass)
+    
+    Float_t Beta_x_Gamma    = fTrackMomentum.P()/mass;//  pc/mc2
+    Float_t SigmaEffect_10degrees;
+    Float_t SigmaEffect_thetadegrees;
+    Float_t ELossParticle_ELossMip;
+    Float_t YAngleEffect=0.;
+    Float_t theta_wires      =  TMath::Abs( TMath::ASin( TMath::Sin(theta) * TMath::Sin(phi) ) );
+    
+    if ( (Beta_x_Gamma >3.2)   &&  (theta_wires*kRaddeg<=10) ) {
+      Beta_x_Gamma=TMath::Log(Beta_x_Gamma);
+      ELossParticle_ELossMip = fElossRatio->Eval(Beta_x_Gamma);
+      // 10 degrees is a reference for a model (arbitrary)
+      SigmaEffect_10degrees=fAngleEffect10->Eval(ELossParticle_ELossMip);// in micrometers
+      // Angle with respect to the wires assuming that chambers are perpendicular to the z axis.
+      SigmaEffect_thetadegrees =  SigmaEffect_10degrees/fAngleEffectNorma->Eval(theta_wires*kRaddeg);  // For 5mm gap  
+      if ( (iChamber==1)  ||  (iChamber==2) )  
+       SigmaEffect_thetadegrees/=(1.09833e+00+1.70000e-02*theta_wires*kRaddeg); // The gap is different (4mm)
+      YAngleEffect=1.e-04*gRandom->Gaus(0,SigmaEffect_thetadegrees); // Error due to the angle effect in cm
+    }
+    
+    
+    // One hit per chamber
+    GetMUONData()->AddHit(fIshunt, gAlice->GetCurrentTrackNumber(), iChamber, ipart, 
+                         fTrackPosition.X(), fTrackPosition.Y()+YAngleEffect, fTrackPosition.Z(), 0.0, 
+                         fTrackMomentum.P(),theta, phi, fStepSum[idvol], fDestepSum[idvol],
+                         fTrackPosition.X(),fTrackPosition.Y(),fTrackPosition.Z());
+//     if (GetDebug()){
+//       Info("StepManager Exit","Particle exiting from chamber %d",iChamber);
+//       Info("StepManager Exit","StepSum %f eloss geant %g ",fStepSum[idvol],fDestepSum[idvol]);
+//       Info("StepManager Exit","Track Position %f %f %f",fTrackPosition.X(),fTrackPosition.Y(),fTrackPosition.Z()) ;
+//     }
+    fStepSum[idvol]  =0; // Reset for the next event
+    fDestepSum[idvol]=0; // Reset for the next event
+  }
 }
 
 //___________________________________________
index a59b1a6842f55522c2d53ab7929f3c230b170b35..02ca246ae7a84139232b00d0c9605b3f6d76e7cf 100644 (file)
 //  Manager and hits classes for set:MUON version 0    //
 /////////////////////////////////////////////////////////
  
+#include "TLorentzVector.h"
+
 #include "AliMUON.h"
-#include <TF1.h>
+
+class TF1;
 
 class AliMUONv1 : public AliMUON {
 public:
@@ -24,28 +27,29 @@ public:
    virtual Int_t  IsVersion() const {return 1;}
    virtual void   StepManager();
    void StepManagerOld();
-   void StepManagerNew();
-   void StepManagerTest();
-
-
-
    void SetStepManagerVersionOld(Bool_t Opt) 
      { fStepManagerVersionOld = Opt; }
-   void SetStepManagerVersionNew(Bool_t Opt) 
-     { fStepManagerVersionNew = Opt; }
-   void SetStepManagerVersionTest(Bool_t Opt) 
-     { fStepManagerVersionTest = Opt; }
    void SetStepMaxInActiveGas(Float_t StepMax)
      {fStepMaxInActiveGas = StepMax; }
 protected:
-   Int_t*  fStations; //! allow to externally set which station to create
+   Int_t*  fStations;              //! allow to externally set which station to create
    Bool_t  fStepManagerVersionOld; // Version of StepManager, Default is false
-   Bool_t  fStepManagerVersionNew; // Version of StepManager, Default is false
-   Bool_t  fStepManagerVersionTest; // Version of StepManager, Default is false
-   Float_t fStepMaxInActiveGas; // Step mas in active gas default 0.6cm
+   Float_t fStepMaxInActiveGas;    // Step max in active gas default 0.6cm
    virtual Int_t  GetChamberId(Int_t volId) const;
-   
 
+   // StepManager 
+   Float_t *  fStepSum; //!
+   Float_t *  fDestepSum; //!
+   // Momentum of the particle entering in the active gas of chamber
+   TLorentzVector fTrackMomentum;
+   // Position of the particle exiting the active gas of chamber
+   TLorentzVector fTrackPosition;
+   // Ratio of particle mean eloss with respect MIP's
+   TF1 * fElossRatio;
+   // Angle effect in tracking chambers at theta =10 degres as a function of ElossRatio (Khalil BOUDJEMLINE sep 2003 Ph.D Thesis) (in micrometers)
+   TF1 * fAngleEffect10;
+   // Angle effect: Normalisation form theta=10 degres to theta between 0 and 10 (Khalil BOUDJEMLINE sep 2003 Ph.D Thesis) 
+   TF1 * fAngleEffectNorma;
 private:
    ClassDef(AliMUONv1,1)  // MUON Detector class Version 1
 
index 1a8642b5c8ea855d4df8a8ea81c1fec9f1547859..81a5d7852fb06d5d27939f1633c86d26a30c7ff0 100644 (file)
@@ -6,7 +6,7 @@ void Config(char directory[100]="", char option[6]="box")
 {
   //
   // Config file for MUON test
-  // Gines MARITNEZ, Subatech, mai 2003, august 2003
+  // Gines MARTINEZ, Subatech, mai 2003, august 2003
   // 
 
   //=====================================================================
@@ -87,8 +87,8 @@ void Config(char directory[100]="", char option[6]="box")
   if (!strcmp(option,"box")) {
     AliGenBox * gener = new AliGenBox(1);
     gener->SetMomentumRange(20.,20.1);
-    gener->SetPhiRange(45., 45.01);         
-    gener->SetThetaRange(4.000,4.001);
+    gener->SetPhiRange(-180, 180.01);         
+    gener->SetThetaRange(171.000,178.001);
     gener->SetPart(13);           // Muons
     gener->SetOrigin(0.,0., 0.);  //vertex position
     gener->SetSigma(0.0, 0.0, 0.0);         //Sigma in (X,Y,Z) (cm) on IP position
@@ -112,9 +112,8 @@ void Config(char directory[100]="", char option[6]="box")
     gener->SetMomentumRange(0,999);
     gener->SetPtRange(0,100.);
     gener->SetPhiRange(-180, 180);
-    gener->SetYRange(2.5,4);
     gener->SetCutOnChild(1);
-    gener->SetChildThetaRange(2.0,9);
+    gener->SetChildThetaRange(171.0,178.0);
     gener->SetOrigin(0,0,0);          //vertex position    gener->SetSigma(0,0,0);           //Sigma in (X,Y,Z) (cm) on IP position
     gener->SetForceDecay(kDiMuon);
     gener->SetTrackingFlag(1);
@@ -123,21 +122,18 @@ void Config(char directory[100]="", char option[6]="box")
      
   //============================================================= 
   // Field (L3 0.4 T)
-   AliMagFMaps* field = new AliMagFMaps("Maps","Maps", 2, 1., 10., AliMagFMaps::k4kG);
-   gAlice->SetField(field);
-  //gAlice->SetField(2,1) ;  //(-999,2);
-
+  AliMagFMaps* field = new AliMagFMaps("Maps","Maps", 2, 1., 10., AliMagFMaps::k4kG);
+  gAlice->SetField(field);
 
   //=================== Alice BODY parameters =============================
   AliBODY *BODY = new AliBODY("BODY","Alice envelop");
 
 
   //=================== MUON Subsystem ===========================
-  AliMUONv1 *MUON  = new AliMUONv1("MUON","default");
-  
+   cout << ">>> Config_MUON_test.C: Creating AliMUONv1 ..."<<endl;
+   AliMUONv1 *MUON  = new AliMUONv1("MUON","default"); 
 }
 
-
 Float_t EtaToTheta(Float_t arg){
   return (180./TMath::Pi())*2.*atan(exp(-arg));
 }
index 52efe77fd0301a4b89436ee19a5e532d60c32cfe..c3d790f95a27bac315e049a1284f483f6771e8ac 100644 (file)
@@ -91,6 +91,16 @@ void MUONrawclusters (char* filename="galice.root", Int_t evNumber1=0,Int_t evNu
 
   MUONLoader->LoadDigits("READ");
   MUONLoader->LoadRecPoints("UPDATE");
+
+  // Testing if RawClusterisation has already been done
+  RunLoader->GetEvent(0);
+  if (MUONLoader->TreeR()) {
+    if (muondata->IsRawClusterBranchesInTree()) {
+      MUONLoader->UnloadRecPoints();
+      MUONLoader->LoadRecPoints("RECREATE");
+      printf("Recreating RecPoints files\n");
+    }
+  }
 //
 //   Loop over events              
   //