Fork me on GitHub

Changes in / [740d430:4a0d9d5] in git


Ignore:
Files:
11 added
14 edited

Legend:

Unmodified
Added
Removed

  • Makefile

    r740d430 r4a0d9d5  
    316316        modules/EnergySmearing.h \
    317317        modules/MomentumSmearing.h \
     318        modules/TrackSmearing.h \
    318319        modules/ImpactParameterSmearing.h \
    319320        modules/TimeSmearing.h \
    … …  
    338339        modules/StatusPidFilter.h \
    339340        modules/PdgCodeFilter.h \
     341        modules/BeamSpotFilter.h \
    340342        modules/Cloner.h \
    341343        modules/Weighter.h \
    … …  
    343345        modules/JetFlavorAssociation.h \
    344346        modules/JetFakeParticle.h \
     347        modules/VertexSorter.h \
     348        modules/VertexFinder.h \
    345349        modules/ExampleModule.h
    346350ModulesDict$(PcmSuf): \
    … …  
    546550        classes/DelphesFactory.h \
    547551        classes/DelphesFormula.h
     552tmp/modules/BeamSpotFilter.$(ObjSuf): \
     553        modules/BeamSpotFilter.$(SrcSuf) \
     554        modules/BeamSpotFilter.h \
     555        classes/DelphesClasses.h \
     556        classes/DelphesFactory.h \
     557        classes/DelphesFormula.h \
     558        external/ExRootAnalysis/ExRootResult.h \
     559        external/ExRootAnalysis/ExRootFilter.h \
     560        external/ExRootAnalysis/ExRootClassifier.h
    548561tmp/modules/Calorimeter.$(ObjSuf): \
    549562        modules/Calorimeter.$(SrcSuf) \
    … …  
    836849        external/ExRootAnalysis/ExRootFilter.h \
    837850        external/ExRootAnalysis/ExRootClassifier.h
     851tmp/modules/TrackSmearing.$(ObjSuf): \
     852        modules/TrackSmearing.$(SrcSuf) \
     853        modules/TrackSmearing.h \
     854        classes/DelphesClasses.h \
     855        classes/DelphesFactory.h \
     856        classes/DelphesFormula.h \
     857        external/ExRootAnalysis/ExRootResult.h \
     858        external/ExRootAnalysis/ExRootFilter.h \
     859        external/ExRootAnalysis/ExRootClassifier.h
    838860tmp/modules/TreeWriter.$(ObjSuf): \
    839861        modules/TreeWriter.$(SrcSuf) \
    … …  
    852874        classes/DelphesFactory.h \
    853875        classes/DelphesFormula.h \
     876        external/ExRootAnalysis/ExRootResult.h \
     877        external/ExRootAnalysis/ExRootFilter.h \
     878        external/ExRootAnalysis/ExRootClassifier.h
     879tmp/modules/VertexFinder.$(ObjSuf): \
     880        modules/VertexFinder.$(SrcSuf) \
     881        modules/VertexFinder.h \
     882        classes/DelphesClasses.h \
     883        classes/DelphesFactory.h \
     884        classes/DelphesFormula.h \
     885        classes/DelphesPileUpReader.h \
     886        external/ExRootAnalysis/ExRootResult.h \
     887        external/ExRootAnalysis/ExRootFilter.h \
     888        external/ExRootAnalysis/ExRootClassifier.h
     889tmp/modules/VertexSorter.$(ObjSuf): \
     890        modules/VertexSorter.$(SrcSuf) \
     891        modules/VertexSorter.h \
     892        classes/DelphesClasses.h \
     893        classes/DelphesFactory.h \
     894        classes/DelphesFormula.h \
     895        classes/DelphesPileUpReader.h \
    854896        external/ExRootAnalysis/ExRootResult.h \
    855897        external/ExRootAnalysis/ExRootFilter.h \
    … …  
    915957        tmp/modules/AngularSmearing.$(ObjSuf) \
    916958        tmp/modules/BTagging.$(ObjSuf) \
     959        tmp/modules/BeamSpotFilter.$(ObjSuf) \
    917960        tmp/modules/Calorimeter.$(ObjSuf) \
    918961        tmp/modules/Cloner.$(ObjSuf) \
    … …  
    946989        tmp/modules/TrackCountingTauTagging.$(ObjSuf) \
    947990        tmp/modules/TrackPileUpSubtractor.$(ObjSuf) \
     991        tmp/modules/TrackSmearing.$(ObjSuf) \
    948992        tmp/modules/TreeWriter.$(ObjSuf) \
    949993        tmp/modules/UniqueObjectFinder.$(ObjSuf) \
     994        tmp/modules/VertexFinder.$(ObjSuf) \
     995        tmp/modules/VertexSorter.$(ObjSuf) \
    950996        tmp/modules/Weighter.$(ObjSuf)
    951997
    … …  
    15521598        tmp/external/tcl/tclVar.$(ObjSuf)
    15531599
     1600modules/TrackSmearing.h: \
     1601        classes/DelphesModule.h
     1602        @touch $@
     1603
    15541604external/fastjet/ClusterSequence.hh: \
    15551605        external/fastjet/PseudoJet.hh \
    … …  
    18091859        @touch $@
    18101860
     1861modules/VertexSorter.h: \
     1862        classes/DelphesModule.h \
     1863        classes/DelphesClasses.h
     1864        @touch $@
     1865
    18111866modules/Delphes.h: \
    18121867        classes/DelphesModule.h
     1868        @touch $@
     1869
     1870modules/VertexFinder.h: \
     1871        classes/DelphesModule.h \
     1872        classes/DelphesClasses.h
    18131873        @touch $@
    18141874
    … …  
    19772037
    19782038modules/FastJetFinder.h: \
     2039        classes/DelphesModule.h
     2040        @touch $@
     2041
     2042modules/BeamSpotFilter.h: \
    19792043        classes/DelphesModule.h
    19802044        @touch $@

  • classes/DelphesClasses.cc

    r740d430 r4a0d9d5  
    4141CompBase *Tower::fgCompare = CompE<Tower>::Instance();
    4242CompBase *HectorHit::fgCompare = CompE<HectorHit>::Instance();
     43CompBase *Vertex::fgCompare = CompSumPT2<Vertex>::Instance();
    4344CompBase *Candidate::fgCompare = CompMomentumPt<Candidate>::Instance();
    4445
    … …  
    121122  Charge(0), Mass(0.0),
    122123  IsPU(0), IsRecoPU(0), IsConstituent(0), IsFromConversion(0),
     124  ClusterIndex(-1), ClusterNDF(0), ClusterSigma(0), SumPT2(0), BTVSumPT2(0), GenDeltaZ(0), GenSumPT2(0),
    123125  Flavor(0), FlavorAlgo(0), FlavorPhys(0),
    124126  BTag(0), BTagAlgo(0), BTagPhys(0),
    … …  
    127129  Momentum(0.0, 0.0, 0.0, 0.0),
    128130  Position(0.0, 0.0, 0.0, 0.0),
     131  PositionError(0.0, 0.0, 0.0, 0.0),
     132  InitialPosition(0.0, 0.0, 0.0, 0.0),
    129133  Area(0.0, 0.0, 0.0, 0.0),
    130   Dxy(0), SDxy(0), Xd(0), Yd(0), Zd(0),
     134  L(0),
     135  D0(0), ErrorD0(0),
     136  DZ(0), ErrorDZ(0),
     137  P(0),  ErrorP(0),
     138  PT(0), ErrorPT(0),
     139  CtgTheta(0), ErrorCtgTheta(0),
     140  Phi(0), ErrorPhi(0), 
     141  Xd(0), Yd(0), Zd(0),
    131142  TrackResolution(0),
    132143  NCharged(0),
    … …  
    245256  object.IsConstituent = IsConstituent;
    246257  object.IsFromConversion = IsFromConversion;
     258  object.ClusterIndex = ClusterIndex;
     259  object.ClusterNDF = ClusterNDF;
     260  object.ClusterSigma = ClusterSigma;
     261  object.SumPT2 = SumPT2;
     262  object.BTVSumPT2 = BTVSumPT2;
     263  object.GenDeltaZ = GenDeltaZ;
     264  object.GenSumPT2 = GenSumPT2;
    247265  object.Flavor = Flavor;
    248266  object.FlavorAlgo = FlavorAlgo;
    … …  
    262280  object.Momentum = Momentum;
    263281  object.Position = Position;
     282  object.InitialPosition = InitialPosition;
     283  object.PositionError = PositionError;
    264284  object.Area = Area;
    265   object.Dxy = Dxy;
    266   object.SDxy = SDxy;
     285  object.L = L;
     286  object.ErrorT = ErrorT;
     287  object.D0 = D0;
     288  object.ErrorD0 = ErrorD0;
     289  object.DZ = DZ;
     290  object.ErrorDZ = ErrorDZ;
     291  object.P = P;
     292  object.ErrorP = ErrorP;
     293  object.PT = PT;
     294  object.ErrorPT = ErrorPT;
     295  object.CtgTheta = CtgTheta ;
     296  object.ErrorCtgTheta = ErrorCtgTheta;
     297  object.Phi = Phi;
     298  object.ErrorPhi = ErrorPhi; 
    267299  object.Xd = Xd;
    268300  object.Yd = Yd;
    … …  
    282314  object.SumPtChargedPU = SumPtChargedPU;
    283315  object.SumPt = SumPt;
    284 
     316  object.ClusterIndex = ClusterIndex;
     317  object.ClusterNDF = ClusterNDF;
     318  object.ClusterSigma = ClusterSigma;
     319  object.SumPT2 = SumPT2;
     320 
    285321  object.FracPt[0] = FracPt[0];
    286322  object.FracPt[1] = FracPt[1];
    … …  
    363399  Momentum.SetXYZT(0.0, 0.0, 0.0, 0.0);
    364400  Position.SetXYZT(0.0, 0.0, 0.0, 0.0);
     401  InitialPosition.SetXYZT(0.0, 0.0, 0.0, 0.0);
    365402  Area.SetXYZT(0.0, 0.0, 0.0, 0.0);
    366   Dxy = 0.0;
    367   SDxy = 0.0;
     403  L = 0.0;
     404  ErrorT = 0.0;
     405  D0 = 0.0; 
     406  ErrorD0 = 0.0;
     407  DZ = 0.0;
     408  ErrorDZ = 0.0;
     409  P =0.0;
     410  ErrorP =0.0;
     411  PT = 0.0;
     412  ErrorPT = 0.0;
     413  CtgTheta = 0.0;
     414  ErrorCtgTheta = 0.0;
     415  Phi = 0.0;
     416  ErrorPhi = 0.0;
    368417  Xd = 0.0;
    369418  Yd = 0.0;
    … …  
    387436  SumPt = -999;
    388437
     438  ClusterIndex = -1;
     439  ClusterNDF = -99;
     440  ClusterSigma = 0.0;
     441  SumPT2 = 0.0;
     442  BTVSumPT2 = 0.0;
     443  GenDeltaZ = 0.0;
     444  GenSumPT2 = 0.0;
     445 
    389446  FracPt[0] = 0.0;
    390447  FracPt[1] = 0.0;

  • classes/DelphesClasses.h

    r740d430 r4a0d9d5  
    147147  Float_t Pz; // particle momentum vector (z component) | hepevt.phep[number][2]
    148148
    149   Float_t PT; // particle transverse momentum
     149  Float_t D0;
     150  Float_t DZ;
     151  Float_t P;
     152  Float_t PT;
     153  Float_t CtgTheta;
     154  Float_t Phi;
    150155  Float_t Eta; // particle pseudorapidity
    151   Float_t Phi; // particle azimuthal angle
    152 
    153156  Float_t Rapidity; // particle rapidity
    154157
    … …  
    168171//---------------------------------------------------------------------------
    169172
    170 class Vertex: public TObject
     173class Vertex: public SortableObject
    171174{
    172175public:
    … …  
    175178  Float_t Y; // vertex position (y component)
    176179  Float_t Z; // vertex position (z component)
    177 
    178   ClassDef(Vertex, 1)
     180 
     181  Double_t ErrorX;
     182  Double_t ErrorY;
     183  Double_t ErrorZ;
     184
     185  Int_t Index;
     186  Int_t NDF;
     187  Double_t Sigma;
     188  Double_t SumPT2;
     189  Double_t BTVSumPT2;
     190  Double_t GenDeltaZ;
     191  Double_t GenSumPT2;
     192
     193  static CompBase *fgCompare; //!
     194  const CompBase *GetCompare() const { return fgCompare; }
     195
     196  ClassDef(Vertex, 3)
    179197};
    180198
    … …  
    397415  Int_t Charge; // track charge
    398416
    399   Float_t PT; // track transverse momentum
    400 
    401417  Float_t Eta; // track pseudorapidity
    402   Float_t Phi; // track azimuthal angle
    403 
     418 
    404419  Float_t EtaOuter; // track pseudorapidity at the tracker edge
    405420  Float_t PhiOuter; // track azimuthal angle at the tracker edge
    … …  
    415430  Float_t TOuter; // track position (z component) at the tracker edge
    416431
    417   Float_t Dxy;     // track signed transverse impact parameter
    418   Float_t SDxy;    // signed error on the track signed transverse impact parameter
     432  Float_t L; // track path length
     433  Float_t ErrorT; // error on the time measurement
     434     
     435  Float_t D0;     // track signed transverse impact parameter
     436  Float_t ErrorD0;    // signed error on the track signed transverse impact parameter
     437 
     438  Float_t DZ; // track transverse momentum
     439  Float_t ErrorDZ; // track transverse momentum error
     440 
     441  Float_t P; // track transverse momentum
     442  Float_t ErrorP; // track transverse momentum error
     443 
     444  Float_t PT; // track transverse momentum
     445  Float_t ErrorPT; // track transverse momentum error
     446 
     447  Float_t CtgTheta; // track transverse momentum
     448  Float_t ErrorCtgTheta; // track transverse momentum error
     449 
     450  Float_t Phi; // track azimuthal angle
     451  Float_t ErrorPhi; // track azimuthal angle
     452 
    419453  Float_t Xd;      // X coordinate of point of closest approach to vertex
    420454  Float_t Yd;      // Y coordinate of point of closest approach to vertex
    … …  
    423457  TRef Particle; // reference to generated particle
    424458
     459  Int_t VertexIndex; // reference to vertex
     460 
    425461  static CompBase *fgCompare; //!
    426462  const CompBase *GetCompare() const { return fgCompare; }
    … …  
    526562  Float_t DeltaPhi;
    527563
    528   TLorentzVector Momentum, Position, Area;
    529 
    530   Float_t Dxy;
    531   Float_t SDxy;
     564  TLorentzVector Momentum, Position, InitialPosition, PositionError, Area;
     565
     566  Float_t L; // path length
     567  Float_t ErrorT; // path length
     568  Float_t D0;
     569  Float_t ErrorD0;
     570  Float_t DZ;
     571  Float_t ErrorDZ;
     572  Float_t P;
     573  Float_t ErrorP;
     574  Float_t PT;
     575  Float_t ErrorPT;
     576  Float_t CtgTheta;
     577  Float_t ErrorCtgTheta;
     578  Float_t Phi;
     579  Float_t ErrorPhi;
     580
    532581  Float_t Xd;
    533582  Float_t Yd;
    … …  
    562611  Float_t SumPt;
    563612
     613  // vertex variables
     614 
     615  Int_t ClusterIndex;
     616  Int_t ClusterNDF;
     617  Double_t ClusterSigma;
     618  Double_t SumPT2;
     619  Double_t BTVSumPT2;
     620  Double_t GenDeltaZ;
     621  Double_t GenSumPT2;
     622
    564623  // N-subjettiness variables
    565624
    … …  
    595654  void SetFactory(DelphesFactory *factory) { fFactory = factory; }
    596655
    597   ClassDef(Candidate, 4)
     656  ClassDef(Candidate, 5)
    598657};
    599658

  • classes/SortableObject.h

    r740d430 r4a0d9d5  
    156156      return -1;
    157157    else if(t1->ET < t2->ET)
     158      return 1;
     159    else
     160      return 0;
     161  }
     162};
     163
     164//---------------------------------------------------------------------------
     165
     166template <typename T>
     167class CompSumPT2: public CompBase
     168{
     169  CompSumPT2() {}
     170public:
     171  static CompSumPT2 *Instance()
     172  {
     173    static CompSumPT2 single;
     174    return &single;
     175  }
     176
     177  Int_t Compare(const TObject *obj1, const TObject *obj2) const
     178  {
     179    const T *t1 = static_cast<const T*>(obj1);
     180    const T *t2 = static_cast<const T*>(obj2);
     181    if(t1->SumPT2 > t2->SumPT2)
     182      return -1;
     183    else if(t1->SumPT2 < t2->SumPT2)
    158184      return 1;
    159185    else

  • doc/genMakefile.tcl

    r740d430 r4a0d9d5  
    260260executableDeps {converters/*.cpp} {examples/*.cpp}
    261261
    262 executableDeps {readers/DelphesHepMC.cpp} {readers/DelphesLHEF.cpp} {readers/DelphesSTDHEP.cpp}
     262executableDeps {readers/DelphesHepMC.cpp} {readers/DelphesLHEF.cpp} {readers/DelphesSTDHEP.cpp} {readers/DelphesROOT.cpp}
    263263
    264264puts {ifeq ($(HAS_CMSSW),true)}

  • modules/ImpactParameterSmearing.cc

    r740d430 r4a0d9d5  
    9696{
    9797  Candidate *candidate, *particle, *mother;
    98   Double_t xd, yd, zd, dxy, sx, sy, sz, ddxy;
     98  Double_t xd, yd, zd, d0, sx, sy, sz, dd0;
    9999  Double_t pt, eta, px, py, phi, e;
    100100
    … …  
    103103  {
    104104
    105     // take momentum before smearing (otherwise apply double smearing on dxy)
     105    // take momentum before smearing (otherwise apply double smearing on d0)
    106106    particle = static_cast<Candidate*>(candidate->GetCandidates()->At(0));
    107107
    … …  
    131131
    132132    // calculate impact parameter (after-smearing)
    133     dxy = (xd*py - yd*px)/pt;
     133    d0 = (xd*py - yd*px)/pt;
    134134
    135     ddxy = gRandom->Gaus(0.0, fFormula->Eval(pt, eta, phi, e));
     135    dd0 = gRandom->Gaus(0.0, fFormula->Eval(pt, eta, phi, e));
    136136
    137137    // fill smeared values in candidate
    … …  
    143143    candidate->Zd = zd;
    144144
    145     candidate->Dxy = dxy;
    146     candidate->SDxy = ddxy;
     145    candidate->D0 = d0;
     146    candidate->ErrorD0 = dd0;
    147147
    148148    candidate->AddCandidate(mother);

  • modules/ModulesLinkDef.h

    r740d430 r4a0d9d5  
    3535#include "modules/EnergySmearing.h"
    3636#include "modules/MomentumSmearing.h"
     37#include "modules/TrackSmearing.h"
    3738#include "modules/ImpactParameterSmearing.h"
    3839#include "modules/TimeSmearing.h"
    … …  
    5758#include "modules/StatusPidFilter.h"
    5859#include "modules/PdgCodeFilter.h"
     60#include "modules/BeamSpotFilter.h"
    5961#include "modules/Cloner.h"
    6062#include "modules/Weighter.h"
    … …  
    6264#include "modules/JetFlavorAssociation.h"
    6365#include "modules/JetFakeParticle.h"
     66#include "modules/VertexSorter.h"
     67#include "modules/VertexFinder.h"
     68#include "modules/VertexFinder4D.h"
    6469#include "modules/ExampleModule.h"
    6570
    … …  
    7984#pragma link C++ class EnergySmearing+;
    8085#pragma link C++ class MomentumSmearing+;
     86#pragma link C++ class TrackSmearing+;
    8187#pragma link C++ class ImpactParameterSmearing+;
    8288#pragma link C++ class TimeSmearing+;
    … …  
    101107#pragma link C++ class StatusPidFilter+;
    102108#pragma link C++ class PdgCodeFilter+;
     109#pragma link C++ class BeamSpotFilter+;
    103110#pragma link C++ class Cloner+;
    104111#pragma link C++ class Weighter+;
    … …  
    106113#pragma link C++ class JetFlavorAssociation+;
    107114#pragma link C++ class JetFakeParticle+;
     115#pragma link C++ class VertexSorter+;
     116#pragma link C++ class VertexFinder+;
     117#pragma link C++ class VertexFinder4D+;
    108118#pragma link C++ class ExampleModule+;
    109119

  • modules/ParticlePropagator.cc

    r740d430 r4a0d9d5  
    6767}
    6868
     69
    6970//------------------------------------------------------------------------------
    7071
    … …  
    9192  fItInputArray = fInputArray->MakeIterator();
    9293
     94  // import beamspot
     95  try
     96  {
     97    fBeamSpotInputArray = ImportArray(GetString("BeamSpotInputArray", "BeamSpotFilter/beamSpotParticle"));
     98  }
     99  catch(runtime_error &e)
     100  {
     101    fBeamSpotInputArray = 0;
     102  } 
    93103  // create output arrays
    94104
    … …  
    111121{
    112122  Candidate *candidate, *mother;
    113   TLorentzVector candidatePosition, candidateMomentum;
     123  TLorentzVector candidatePosition, candidateMomentum, beamSpotPosition;
    114124  Double_t px, py, pz, pt, pt2, e, q;
    115125  Double_t x, y, z, t, r, phi;
    … …  
    120130  Double_t tmp, discr, discr2;
    121131  Double_t delta, gammam, omega, asinrho;
    122   Double_t rcu, rc2, dxy, xd, yd, zd;
    123 
     132  Double_t rcu, rc2, xd, yd, zd;
     133  Double_t l, d0, dz, p, ctgTheta, phip, etap, alpha;
     134  Double_t bsx, bsy, bsz;
     135         
    124136  const Double_t c_light = 2.99792458E8;
    125 
     137 
     138  if (!fBeamSpotInputArray || fBeamSpotInputArray->GetSize () == 0)
     139    beamSpotPosition.SetXYZT(0.0, 0.0, 0.0, 0.0);
     140  else
     141  {
     142    Candidate &beamSpotCandidate = *((Candidate *) fBeamSpotInputArray->At(0));
     143    beamSpotPosition = beamSpotCandidate.Position;
     144  }
     145 
    126146  fItInputArray->Reset();
    127147  while((candidate = static_cast<Candidate*>(fItInputArray->Next())))
    … …  
    132152    y = candidatePosition.Y()*1.0E-3;
    133153    z = candidatePosition.Z()*1.0E-3;
     154   
     155    bsx = beamSpotPosition.X()*1.0E-3;
     156    bsy = beamSpotPosition.Y()*1.0E-3;
     157    bsz = beamSpotPosition.Z()*1.0E-3;
     158   
    134159    q = candidate->Charge;
    135160
    … …  
    181206      y_t = y + py*t;
    182207      z_t = z + pz*t;
     208     
     209      l = TMath::Sqrt( (x_t - x)*(x_t - x) + (y_t - y)*(y_t - y) + (z_t - z)*(z_t - z));
    183210
    184211      mother = candidate;
    185212      candidate = static_cast<Candidate*>(candidate->Clone());
    186213
     214      candidate->InitialPosition = candidatePosition;
    187215      candidate->Position.SetXYZT(x_t*1.0E3, y_t*1.0E3, z_t*1.0E3, candidatePosition.T() + t*e*1.0E3);
    188 
     216      candidate->L = l*1.0E3;
     217   
    189218      candidate->Momentum = candidateMomentum;
    190219      candidate->AddCandidate(mother);
    … …  
    238267      yd = (rc2 > 0.0) ? yd / rc2 : -999;
    239268      zd = z + (TMath::Sqrt(xd*xd + yd*yd) - TMath::Sqrt(x*x + y*y))*pz/pt;
    240 
    241       // calculate impact paramater
    242       dxy = (xd*py - yd*px)/pt;
    243 
     269     
     270      // use perigee momentum rather than original particle
     271      // momentum, since the orignal particle momentum isn't known
     272     
     273      px = TMath::Sign(1.0,r) * pt * (-y_c / r_c);
     274      py = TMath::Sign(1.0,r) * pt * (x_c / r_c);
     275      etap = candidateMomentum.Eta();
     276      phip = TMath::ATan2(py, px);
     277       
     278      candidateMomentum.SetPtEtaPhiE(pt, etap, phip, candidateMomentum.E());
     279     
     280      // calculate additional track parameters (correct for beamspot position)
     281       
     282      d0        = (  (x - bsx) * py - (y - bsy) * px) / pt;
     283      dz        = z - ((x - bsx) * px + (y - bsy) * py) / pt * (pz / pt);
     284      p         = candidateMomentum.P();
     285      ctgTheta  = 1.0 / TMath::Tan (candidateMomentum.Theta ());
     286         
    244287      // 3. time evaluation t = TMath::Min(t_r, t_z)
    245288      //    t_r : time to exit from the sides
    … …  
    287330      r_t = TMath::Hypot(x_t, y_t);
    288331
     332     
     333      // compute path length for an helix
     334     
     335      alpha = pz*1.0E9 / c_light / gammam;
     336      l = t * TMath::Sqrt(alpha*alpha + r*r*omega*omega);
     337         
    289338      if(r_t > 0.0)
    290339      {
     340       
     341        // store these variables before cloning
     342           
     343        candidate->D0 = d0*1.0E3;
     344        candidate->DZ = dz*1.0E3;
     345        candidate->P  = p;
     346        candidate->PT = pt;
     347        candidate->CtgTheta = ctgTheta;
     348        candidate->Phi = phip;
     349     
    291350        mother = candidate;
    292351        candidate = static_cast<Candidate*>(candidate->Clone());
    293352
     353        candidate->InitialPosition = candidatePosition;
    294354        candidate->Position.SetXYZT(x_t*1.0E3, y_t*1.0E3, z_t*1.0E3, candidatePosition.T() + t*c_light*1.0E3);
    295355
    296356        candidate->Momentum = candidateMomentum;
    297         candidate->Dxy = dxy*1.0E3;
    298         candidate->Xd = xd*1.0E3;
     357       
     358            candidate->L  =  l*1.0E3;
     359       
     360            candidate->Xd = xd*1.0E3;
    299361        candidate->Yd = yd*1.0E3;
    300362        candidate->Zd = zd*1.0E3;

  • modules/ParticlePropagator.h

    r740d430 r4a0d9d5  
    3535class TClonesArray;
    3636class TIterator;
     37class TLorentzVector;
    3738
    3839class ParticlePropagator: public DelphesModule
    … …  
    5556
    5657  const TObjArray *fInputArray; //!
     58  const TObjArray *fBeamSpotInputArray; //!
    5759
    5860  TObjArray *fOutputArray; //!

  • modules/PileUpMerger.cc

    r740d430 r4a0d9d5  
    115115  TDatabasePDG *pdg = TDatabasePDG::Instance();
    116116  TParticlePDG *pdgParticle;
    117   Int_t pid;
     117  Int_t pid, nch, nvtx = -1;
    118118  Float_t x, y, z, t, vx, vy;
    119   Float_t px, py, pz, e;
    120   Double_t dz, dphi, dt;
     119  Float_t px, py, pz, e, pt;
     120  Double_t dz, dphi, dt, sumpt2, dz0, dt0;
    121121  Int_t numberOfEvents, event, numberOfParticles;
    122122  Long64_t allEntries, entry;
    … …  
    132132  fFunction->GetRandom2(dz, dt);
    133133
     134  dz0 = -1.0e6;
     135  dt0 = -1.0e6;
     136 
    134137  dt *= c_light*1.0E3; // necessary in order to make t in mm/c
    135138  dz *= 1.0E3; // necessary in order to make z in mm
     139 
     140  //cout<<dz<<","<<dt<<endl;
     141 
    136142  vx = 0.0;
    137143  vy = 0.0;
     144 
    138145  numberOfParticles = fInputArray->GetEntriesFast();
     146  nch = 0;
     147  sumpt2 = 0.0;
     148   
    139149  while((candidate = static_cast<Candidate*>(fItInputArray->Next())))
    140150  {
    … …  
    143153    z = candidate->Position.Z();
    144154    t = candidate->Position.T();
    145     candidate->Position.SetZ(z + dz);
    146     candidate->Position.SetT(t + dt);
     155    pt = candidate->Momentum.Pt();
     156       
     157    // take postion and time from first stable particle
     158    if (dz0 < -999999.0)
     159      dz0 = z;
     160    if (dt0 < -999999.0)
     161      dt0 = t;
     162
     163    // cancel any possible offset in position and time the input file
     164    candidate->Position.SetZ(z - dz0 + dz);
     165    candidate->Position.SetT(t - dt0 + dt);
     166   
    147167    fParticleOutputArray->Add(candidate);
    148   }
    149 
     168 
     169    if(TMath::Abs(candidate->Charge) >  1.0E-9)
     170    {
     171      nch++;   
     172      sumpt2 += pt*pt;
     173    } 
     174  }
     175 
    150176  if(numberOfParticles > 0)
    151177  {
    152     vx /= numberOfParticles;
    153     vy /= numberOfParticles;
    154   }
    155 
     178    vx /= sumpt2;
     179    vy /= sumpt2;
     180  }
     181 
     182  nvtx++;
    156183  factory = GetFactory();
    157184
    158185  vertex = factory->NewCandidate();
    159186  vertex->Position.SetXYZT(vx, vy, dz, dt);
     187  vertex->ClusterIndex = nvtx;
     188  vertex->ClusterNDF = nch;
     189  vertex->SumPT2 = sumpt2;
     190  vertex->GenSumPT2 = sumpt2;
     191 
    160192  fVertexOutputArray->Add(vertex);
    161193
    … …  
    170202      numberOfEvents = gRandom->Integer(2*fMeanPileUp + 1);
    171203      break;
     204    case 2:
     205      numberOfEvents = fMeanPileUp;
     206      break;   
    172207    default:
    173208      numberOfEvents = gRandom->Poisson(fMeanPileUp);
    … …  
    198233    vx = 0.0;
    199234    vy = 0.0;
     235   
    200236    numberOfParticles = 0;
     237    sumpt2 = 0.0;
     238   
    201239    while(fReader->ReadParticle(pid, x, y, z, t, px, py, pz, e))
    202240    {
    … …  
    224262      vx += candidate->Position.X();
    225263      vy += candidate->Position.Y();
     264     
    226265      ++numberOfParticles;
    227 
     266      if(TMath::Abs(candidate->Charge) >  1.0E-9)
     267      {
     268        nch++;   
     269        sumpt2 += pt*pt;
     270      }
     271       
    228272      fParticleOutputArray->Add(candidate);
    229273    }
    … …  
    234278      vy /= numberOfParticles;
    235279    }
     280   
     281    nvtx++;
    236282
    237283    vertex = factory->NewCandidate();
    238284    vertex->Position.SetXYZT(vx, vy, dz, dt);
     285   
     286    vertex->ClusterIndex = nvtx;
     287    vertex->ClusterNDF = nch;
     288    vertex->SumPT2 = sumpt2;
     289    vertex->GenSumPT2 = sumpt2;
     290   
    239291    vertex->IsPU = 1;
    240292
    241293    fVertexOutputArray->Add(vertex);
    242   }
    243 }
    244 
    245 //------------------------------------------------------------------------------
     294   
     295  }
     296}
     297
     298//------------------------------------------------------------------------------

  • modules/TimeSmearing.cc

    r740d430 r4a0d9d5  
    9999  while((candidate = static_cast<Candidate*>(fItInputArray->Next())))
    100100  {
    101     const TLorentzVector &candidatePosition = candidate->Position;
     101    const TLorentzVector &candidatePosition = candidate->InitialPosition;
    102102    t = candidatePosition.T()*1.0E-3/c_light;
    103103
    104104    // apply smearing formula
    105105    t = gRandom->Gaus(t, fTimeResolution);
    106 
     106   
    107107    mother = candidate;
    108108    candidate = static_cast<Candidate*>(candidate->Clone());
    109     candidate->Position.SetT(t*1.0E3*c_light);
     109    candidate->InitialPosition.SetT(t*1.0E3*c_light);
     110    candidate->ErrorT = fTimeResolution*1.0E3*c_light;
    110111
    111112    candidate->AddCandidate(mother);

  • modules/TrackCountingBTagging.cc

    r740d430 r4a0d9d5  
    9696
    9797  Double_t jpx, jpy;
    98   Double_t dr, tpx, tpy, tpt;
    99   Double_t xd, yd, dxy, ddxy, ip, sip;
     98  Double_t dr, tpt;
     99  Double_t xd, yd, d0, dd0, ip, sip;
    100100
    101101  Int_t sign;
    … …  
    117117    {
    118118      const TLorentzVector &trkMomentum = track->Momentum;
    119 
     119     
    120120      dr = jetMomentum.DeltaR(trkMomentum);
    121 
    122121      tpt = trkMomentum.Pt();
    123       tpx = trkMomentum.Px();
    124       tpy = trkMomentum.Py();
    125 
    126122      xd = track->Xd;
    127123      yd = track->Yd;
    128       dxy = TMath::Hypot(xd, yd);
    129       ddxy = track->SDxy;
     124      d0 = TMath::Hypot(xd, yd);
     125      dd0 = track->ErrorD0;
    130126
    131127      if(tpt < fPtMin) continue;
    132128      if(dr > fDeltaR) continue;
    133       if(dxy > fIPmax) continue;
     129      if(d0 > fIPmax) continue;
    134130
    135131      sign = (jpx*xd + jpy*yd > 0.0) ? 1 : -1;
    136132
    137       ip = sign*dxy;
    138       sip = ip / TMath::Abs(ddxy);
     133      ip = sign*d0;
     134      sip = ip / TMath::Abs(dd0);
    139135
    140136      if(sip > fSigMin) count++;

  • modules/TreeWriter.cc

    r740d430 r4a0d9d5  
    215215    entry->Pz = momentum.Pz();
    216216
     217    entry->D0            = candidate->D0;
     218    entry->DZ            = candidate->DZ;
     219    entry->P             = candidate->P;
     220    entry->PT            = candidate->PT;
     221    entry->CtgTheta      = candidate->CtgTheta;
     222    entry->Phi           = candidate->Phi;
     223
    217224    entry->Eta = eta;
    218225    entry->Phi = momentum.Phi();
    … …  
    237244
    238245  const Double_t c_light = 2.99792458E8;
    239 
     246 
     247  Double_t x, y, z, t, xError, yError, zError, sigma, sumPT2, btvSumPT2, genDeltaZ, genSumPT2;
     248  UInt_t index, ndf;
     249 
     250  array->Sort();
     251 
    240252  // loop over all vertices
    241253  iterator.Reset();
    242254  while((candidate = static_cast<Candidate*>(iterator.Next())))
    243255  {
    244     const TLorentzVector &position = candidate->Position;
     256   
     257    index = candidate->ClusterIndex;
     258    ndf = candidate->ClusterNDF;
     259    sigma = candidate->ClusterSigma;
     260    sumPT2 = candidate->SumPT2;
     261    btvSumPT2 = candidate->BTVSumPT2;
     262    genDeltaZ = candidate->GenDeltaZ;
     263    genSumPT2 = candidate->GenSumPT2;
     264   
     265    x = candidate->Position.X();
     266    y = candidate->Position.Y();
     267    z = candidate->Position.Z();
     268    t = candidate->Position.T()*1.0E-3/c_light;
     269 
     270    xError = candidate->PositionError.X ();
     271    yError = candidate->PositionError.Y ();
     272    zError = candidate->PositionError.Z ();
    245273
    246274    entry = static_cast<Vertex*>(branch->NewEntry());
    247275
    248     entry->X = position.X();
    249     entry->Y = position.Y();
    250     entry->Z = position.Z();
    251     entry->T = position.T()*1.0E-3/c_light;
    252   }
    253 }
     276    entry->Index = index;
     277    entry->NDF = ndf;
     278    entry->Sigma = sigma;
     279    entry->SumPT2 = sumPT2;
     280    entry->BTVSumPT2 = btvSumPT2;
     281    entry->GenDeltaZ = genDeltaZ;
     282    entry->GenSumPT2 = genSumPT2;
     283   
     284    entry->X = x;
     285    entry->Y = y;
     286    entry->Z = z;
     287    entry->T = t;
     288   
     289    entry->ErrorX = xError;
     290    entry->ErrorY = yError;
     291    entry->ErrorZ = zError;
     292  }
     293}
     294
    254295
    255296//------------------------------------------------------------------------------
    … …  
    291332    entry->ZOuter = position.Z();
    292333    entry->TOuter = position.T()*1.0E-3/c_light;
    293 
    294     entry->Dxy = candidate->Dxy;
    295     entry->SDxy = candidate->SDxy ;
     334 
     335    entry->L = candidate->L;
     336 
     337    entry->D0            = candidate->D0;
     338    entry->ErrorD0       = candidate->ErrorD0;
     339    entry->DZ            = candidate->DZ;
     340    entry->ErrorDZ       = candidate->ErrorDZ;
     341    entry->P             = candidate->P;
     342    entry->ErrorP        = candidate->ErrorP;
     343    entry->PT            = candidate->PT;
     344    entry->ErrorPT       = candidate->ErrorPT;
     345    entry->CtgTheta      = candidate->CtgTheta;
     346    entry->ErrorCtgTheta = candidate->ErrorCtgTheta;
     347    entry->Phi           = candidate->Phi;
     348    entry->ErrorPhi      = candidate->ErrorPhi;
     349   
    296350    entry->Xd = candidate->Xd;
    297351    entry->Yd = candidate->Yd;
    … …  
    307361
    308362    entry->Eta = eta;
    309     entry->Phi = momentum.Phi();
    310     entry->PT = pt;
    311 
     363   
    312364    particle = static_cast<Candidate*>(candidate->GetCandidates()->At(0));
    313365    const TLorentzVector &initialPosition = particle->Position;
    … …  
    319371
    320372    entry->Particle = particle;
     373
     374    entry->VertexIndex = candidate->ClusterIndex;
     375
    321376  }
    322377}

  • readers/DelphesCMSFWLite.cpp

    r740d430 r4a0d9d5  
    240240
    241241    branchEvent = treeWriter->NewBranch("Event", HepMCEvent::Class());
    242     branchRwgt = treeWriter->NewBranch("Rwgt", Weight::Class());
     242    branchRwgt = treeWriter->NewBranch("Weight", Weight::Class());
    243243
    244244    confReader = new ExRootConfReader;
Note: See TracChangeset for help on using the changeset viewer.