FASTSIM/AliFastMuonTrackingRes.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 // Implementation of AliFastResponse for the Muon Spectrometer resolution.
  19 // The response depends on the charge of the muon and
  20 // the background level.
  21 // The class uses the instance of an object of type AliMUONFastTracking to
  22 // obtain the smearing parameters.
  23 // Author: andreas.morsch@cern.ch
  24
  25 #include "AliFastMuonTrackingRes.h"
  26 #include "AliMUONFastTracking.h"
  27 #include <TRandom.h>
  28 #include <TF1.h>
  29
  30 ClassImp(AliFastMuonTrackingRes)
  31
  32
  33 AliFastMuonTrackingRes::AliFastMuonTrackingRes() :
  34     AliFastResponse("Resolution", "Muon Tracking Resolution"),
  35     fBackground(0.),
  36     fCharge(1.),
  37     fFastTracking(0)
  38 {
  39 // Default constructor
  40 }
  41
  42 AliFastMuonTrackingRes::AliFastMuonTrackingRes(const AliFastMuonTrackingRes & res)
  43     :AliFastResponse(res),
  44     fBackground(0.),
  45     fCharge(1.),
  46     fFastTracking(0)
  47 {
  48 // Copy constructor
  49     res.Copy(*this);
  50 }
  51
  52 void AliFastMuonTrackingRes::Init()
  53 {
  54 // Initialisation
  55     fFastTracking = AliMUONFastTracking::Instance();
  56     fFastTracking->Init(fBackground);
  57 }
  58
  59
  60
  61 void AliFastMuonTrackingRes::Evaluate(Float_t   p,  Float_t  theta , Float_t   phi,
  62                                          Float_t& pS,  Float_t& thetaS, Float_t&  phiS)
  63 {
  64 //
  65 // Evaluate Gaussian smearing from given kinematics
  66 //
  67
  68     Double_t meanp    = fFastTracking->MeanP  (p, theta, phi, Int_t(fCharge));
  69     Double_t sigmap   = fFastTracking->SigmaP (p, theta, phi, Int_t(fCharge));
  70     Double_t sigma1p  = fFastTracking->Sigma1P(p, theta, phi, Int_t(fCharge));
  71     Double_t normg2   = fFastTracking->NormG2 (p, theta, phi, Int_t(fCharge));
  72     Double_t meang2   = fFastTracking->MeanG2 (p, theta, phi, Int_t(fCharge));
  73     Double_t sigmag2  = fFastTracking->SigmaG2(p, theta, phi, Int_t(fCharge));
  74
  75     Int_t ip,itheta,iphi;
  76     fFastTracking->GetIpIthetaIphi(p, theta, phi, Int_t(fCharge), ip, itheta, iphi);
  77     TF1* fitp = fFastTracking->GetFitP(ip,itheta,iphi);
  78
  79     Float_t curmeanp   = fitp->GetParameter(0);
  80     Float_t cursigmap  = fitp->GetParameter(1);
  81     Float_t cursigma1p = fitp->GetParameter(2);
  82     Float_t curnormg2  = fitp->GetParameter(3);
  83     Float_t curmeang2  = fitp->GetParameter(4);
  84     Float_t cursigmag2 = fitp->GetParameter(5);
  85     if (curmeanp != meanp || cursigmap != sigmap || cursigma1p != sigma1p ||
  86         curnormg2 != normg2 || curmeang2 != meang2 || cursigmag2 != sigmag2){
  87       printf ("Setting new parameters for ip=%d itheta=%d iphi=%d\n",ip,itheta,iphi);
  88       fitp->SetParameters(meanp,sigmap,sigma1p,normg2,meang2,sigmag2);
  89     }
  90
  91     Double_t meantheta  = fFastTracking->MeanTheta (p, theta, phi, Int_t(fCharge));
  92     Double_t sigmatheta = fFastTracking->SigmaTheta(p, theta, phi, Int_t(fCharge));
  93     Double_t meanphi    = fFastTracking->MeanPhi   (p, theta, phi, Int_t(fCharge));
  94     Double_t sigmaphi   = fFastTracking->SigmaPhi  (p, theta, phi, Int_t(fCharge));
  95
  96     if (sigmatheta<0 || sigmaphi<0)
  97       printf ("bin %d %d %d sigmatheta = %f, sigmaphi = %f\n",
  98               ip,itheta,iphi,sigmatheta,sigmaphi);
  99     // Components different from ip=0 have the RMS bigger than mean
 100     Float_t ptp[3]  =  { 1.219576,-0.354764,-0.690117 };
 101     Float_t ptph[3] =  { 0.977522, 0.016269, 0.023158 };
 102     Float_t pphp[3] =  { 1.303256,-0.464847,-0.869322 };
 103
 104     // Smeared momentum
 105     pS = -1.;
 106     //    Float_t dpmax = 5. + ip * 2.5;
 107     //    Float_t dpmax = 5. + ip * 2;
 108     Float_t dpmax;
 109     if (sigmag2<999.) dpmax = 5. * TMath::Abs(sigmap + sigmag2);
 110     else dpmax = 5. * TMath::Abs(sigmap);
 111     Float_t dp = 100;
 112     while (pS<0 || TMath::Abs(dp)>dpmax) {
 113       pS = p + fitp->GetRandom();
 114       dp = pS - p;
 115     }
 116     // Smeared phi
 117     Float_t sigmaphiold=sigmaphi;
 118     if (ip==0) sigmaphi *= pphp[0] + pphp[1] * dp + pphp[2] * dp*dp;
 119     if (sigmaphi<0.5 * sigmaphiold) sigmaphi = 0.5 * sigmaphiold;
 120     if (sigmaphi>2.  * sigmaphiold) sigmaphi = 2.  * sigmaphiold;
 121     phiS = phi + gRandom->Gaus(meanphi, sigmaphi);
 122     Float_t dphi = phiS - phi;
 123     // Smeared theta
 124     Float_t sigmathetaold=sigmatheta;
 125     if (ip==0) sigmatheta *= ptp[0] + ptp[1] * dp + ptp[2] * dp*dp;
 126     if (ip==0) sigmatheta *= ptph[0] + ptph[1] * dphi + ptph[2] * dphi*dphi;
 127     if (sigmatheta<0.5 * sigmathetaold) sigmatheta = 0.5 * sigmathetaold;
 128     if (sigmatheta>2.  * sigmathetaold) sigmatheta = 2.  * sigmathetaold;
 129     thetaS = theta + gRandom->Gaus(meantheta,sigmatheta);
 130 }
 131
 132 AliFastMuonTrackingRes& AliFastMuonTrackingRes::operator=(const  AliFastMuonTrackingRes& rhs)
 133 {
 134 // Assignment operator
 135     rhs.Copy(*this);
 136     return *this;
 137 }
 138
 139
 140
 141