Fork me on GitHub

Ignore:
File:
1 edited

Legend:

Unmodified
Added
Removed
  • modules/TimeSmearing.cc

    r2b5ff2c r7939c6c  
    5353
    5454TimeSmearing::TimeSmearing() :
    55   fFormula(0), fItInputArray(0)
     55  fItInputArray(0)
    5656{
    57   fFormula = new DelphesFormula;
    5857}
    5958
     
    6160
    6261TimeSmearing::~TimeSmearing()
    63 {
    64   if(fFormula) delete fFormula;
     62
    6563}
    6664
     
    7169  // read resolution formula
    7270
    73   fFormula->Compile(GetString("TimeResolution", "1.0"));
     71  fTimeResolution = GetDouble("TimeResolution", 1.);
    7472
    7573  // import input array
     
    7775  fInputArray = ImportArray(GetString("InputArray", "MuonMomentumSmearing/muons"));
    7876  fItInputArray = fInputArray->MakeIterator();
    79 
    8077  // create output array
    8178
     
    9592{
    9693  Candidate *candidate, *mother;
    97   Double_t ti, tf_smeared, tf, timeResolution;
    98   Double_t pt, eta, phi, e, d0, dz, ctgTheta;
     94  Double_t ti, tf_smeared, tf;
     95  Double_t pt, eta, phi, e, p, l;
     96  Double_t sigma_t, beta_particle;
    9997
    10098
    10199  const Double_t c_light = 2.99792458E8;
     100
     101  cout << " STARTINNNGG ---------->" << endl;
    102102
    103103  fItInputArray->Reset();
     
    114114    pt = candidateMomentum.Pt();
    115115    e = candidateMomentum.E();
    116     d0 = candidate->D0;
    117     dz = candidate->DZ;
    118     ctgTheta = candidate->CtgTheta;
     116    p = candidateMomentum.P();
     117    beta_particle = p/e;
     118    l = candidate->L;
     119    // apply smearing formula
    119120
    120     // apply smearing formula
    121    
    122     timeResolution = fFormula->Eval(pt, eta, phi, e, d0, dz, ctgTheta);
    123     if(fabs(candidate->Position.Eta())<fEtaMax)
     121    if(candidate->Charge != 0)
    124122    {
    125       tf_smeared = tf + timeResolution*gRandom->Gaus(0, 1);
     123      tf_smeared = tf + fTimeResolution*gRandom->Gaus(0, 1);
     124      //mother = candidate;
     125      //candidate = static_cast<Candidate*>(candidate->Clone());  // I am not sure that we need these lines !!!
     126      //candidate->AddCandidate(mother);
     127      candidate->InitialPosition.SetT((100+ti)*1.0E3*c_light);
     128      candidate->Position.SetT(tf_smeared*1.0E3*c_light);
     129      candidate->ErrorT = fTimeResolution*1.0E3*c_light;
     130      fOutputArray->Add(candidate);
    126131    }
    127     else continue;
    128 
    129     // double beta_particle = candidate->Momentum.P()/candidate->Momentum.E();
    130     // ti = tf_smeared - candidate->Ld*1.0E-3/(c_light*beta_particle);
    131 
    132     mother = candidate;
    133     candidate = static_cast<Candidate*>(candidate->Clone());
    134     candidate->AddCandidate(mother);
    135     candidate->InitialPosition.SetT((100+ti)*1.0E3*c_light);
    136     candidate->Position.SetT(tf_smeared*1.0E3*c_light);
    137     candidate->ErrorT = timeResolution*1.0E3*c_light;
    138 
    139     fOutputArray->Add(candidate);
     132    else
     133    {
     134        sigma_t = sqrt(pow(20,2) + pow(150,2)/e);
     135        ti = sigma_t - l*1.0E3/(c_light*beta_particle);   
     136        candidate->InitialPosition.SetT(ti);
     137        candidate->ErrorT = sigma_t*1.0E3*c_light;              // Do we need to sum with 100 like in upside ?
     138        fOutputArray->Add(candidate);
     139    }
     140   
    140141  }
    141142}
Note: See TracChangeset for help on using the changeset viewer.