Fork me on GitHub

Changes in / [5eb063e:a643a7f] in git


Ignore:
Files:
54 added
8 deleted
12 edited

Legend:

Unmodified
Added
Removed

  • Makefile

    r5eb063e ra643a7f  
    8989DISPLAYLIB = libDelphesDisplay.lib
    9090
    91 VERSION = $(shell cat VERSION)
     91VERSION = x.y.z
    9292DISTDIR = Delphes-$(VERSION)
    9393DISTTAR = $(DISTDIR).tar.gz
    … …  
    174174tmp/examples/Example1.$(ObjSuf): \
    175175        examples/Example1.cpp \
    176         classes/DelphesClasses.h \
    177         external/ExRootAnalysis/ExRootResult.h \
    178         external/ExRootAnalysis/ExRootTreeBranch.h \
    179         external/ExRootAnalysis/ExRootTreeReader.h \
    180         external/ExRootAnalysis/ExRootTreeWriter.h \
    181         external/ExRootAnalysis/ExRootUtilities.h
    182 DelphesValidation$(ExeSuf): \
    183         tmp/validation/DelphesValidation.$(ObjSuf)
    184 
    185 tmp/validation/DelphesValidation.$(ObjSuf): \
    186         validation/DelphesValidation.cpp \
    187176        classes/DelphesClasses.h \
    188177        external/ExRootAnalysis/ExRootResult.h \
    … …  
    199188        stdhep2pileup$(ExeSuf) \
    200189        CaloGrid$(ExeSuf) \
    201         Example1$(ExeSuf) \
    202         DelphesValidation$(ExeSuf)
     190        Example1$(ExeSuf)
    203191
    204192EXECUTABLE_OBJ +=  \
    … …  
    210198        tmp/converters/stdhep2pileup.$(ObjSuf) \
    211199        tmp/examples/CaloGrid.$(ObjSuf) \
    212         tmp/examples/Example1.$(ObjSuf) \
    213         tmp/validation/DelphesValidation.$(ObjSuf)
     200        tmp/examples/Example1.$(ObjSuf)
    214201
    215202DelphesHepMC$(ExeSuf): \
    … …  
    24342421        @echo ">> Building $(DISTTAR)"
    24352422        @mkdir -p $(DISTDIR)
    2436         @cp -a AUTHORS CHANGELOG CMakeLists.txt COPYING DelphesEnv.sh LICENSE NOTICE README README_4LHCb VERSION Makefile MinBias.pileup configure cards classes converters display doc examples external modules python readers validation $(DISTDIR)
     2423        @cp -a AUTHORS CHANGELOG CMakeLists.txt COPYING DelphesEnv.sh LICENSE NOTICE README README_4LHCb Makefile MinBias.pileup configure cards classes converters display doc examples external modules python readers validation $(DISTDIR)
    24372424        @find $(DISTDIR) -depth -name .\* -exec rm -rf {} \;
    24382425        @tar -czf $(DISTTAR) $(DISTDIR)

  • classes/DelphesClasses.cc

    r5eb063e ra643a7f  
    2828
    2929#include "classes/DelphesClasses.h"
    30 
    3130#include "classes/DelphesFactory.h"
    3231#include "classes/SortableObject.h"
    … …  
    140139  PositionError(0.0, 0.0, 0.0, 0.0),
    141140  Area(0.0, 0.0, 0.0, 0.0),
     141  TrackCovariance(5),
    142142  L(0),
    143143  D0(0), ErrorD0(0),
    144144  DZ(0), ErrorDZ(0),
    145145  P(0), ErrorP(0),
     146  C(0), ErrorC(0),
    146147  PT(0), ErrorPT(0),
    147148  CtgTheta(0), ErrorCtgTheta(0),
    … …  
    313314  object.P = P;
    314315  object.ErrorP = ErrorP;
     316  object.C = C;
     317  object.ErrorC = ErrorC;
    315318  object.PT = PT;
    316319  object.ErrorPT = ErrorPT;
    … …  
    377380  object.SoftDroppedSubJet1 = SoftDroppedSubJet1;
    378381  object.SoftDroppedSubJet2 = SoftDroppedSubJet2;
    379 
     382  object.TrackCovariance = TrackCovariance;
    380383  object.fFactory = fFactory;
    381384  object.fArray = 0;
    … …  
    434437  InitialPosition.SetXYZT(0.0, 0.0, 0.0, 0.0);
    435438  Area.SetXYZT(0.0, 0.0, 0.0, 0.0);
     439  TrackCovariance.Zero();
    436440  L = 0.0;
    437441  ErrorT = 0.0;
    … …  
    442446  P = 0.0;
    443447  ErrorP = 0.0;
     448  C = 0.0;
     449  ErrorC = 0.0;
    444450  PT = 0.0;
    445451  ErrorPT = 0.0;

  • classes/DelphesClasses.h

    r5eb063e ra643a7f  
    3333
    3434#include "TLorentzVector.h"
     35#include "TMatrixDSym.h"
    3536#include "TObject.h"
    3637#include "TRef.h"
    … …  
    650651
    651652  Float_t L; // path length
     653  Float_t DZ;
     654  Float_t ErrorDZ;
    652655  Float_t ErrorT; // path length
    653656  Float_t D0;
    654657  Float_t ErrorD0;
    655   Float_t DZ;
    656   Float_t ErrorDZ;
     658  Float_t C;
     659  Float_t ErrorC;
    657660  Float_t P;
    658661  Float_t ErrorP;
    … …  
    699702  Float_t SumPt;
    700703
     704  // ACTS compliant 6x6 track covariance (D0, phi, Curvature, dz, ctg(theta))
     705
     706  TMatrixDSym TrackCovariance;
     707
    701708  // vertex variables
    702709

  • doc/genMakefile.tcl

    r5eb063e ra643a7f  
    282282DISPLAYLIB = libDelphesDisplay.lib
    283283
    284 VERSION = $(shell cat VERSION)
     284VERSION = x.y.z
    285285DISTDIR = Delphes-$(VERSION)
    286286DISTTAR = $(DISTDIR).tar.gz
    … …  
    290290}
    291291
    292 executableDeps {converters/*.cpp} {examples/*.cpp} {validation/*.cpp}
     292executableDeps {converters/*.cpp} {examples/*.cpp}
    293293
    294294executableDeps {readers/DelphesHepMC.cpp} {readers/DelphesLHEF.cpp} {readers/DelphesSTDHEP.cpp} {readers/DelphesROOT.cpp}
    … …  
    434434        @echo ">> Building $(DISTTAR)"
    435435        @mkdir -p $(DISTDIR)
    436         @cp -a AUTHORS CHANGELOG CMakeLists.txt COPYING DelphesEnv.sh LICENSE NOTICE README README_4LHCb VERSION Makefile MinBias.pileup configure cards classes converters display doc examples external modules python readers validation $(DISTDIR)
     436        @cp -a AUTHORS CHANGELOG CMakeLists.txt COPYING DelphesEnv.sh LICENSE NOTICE README README_4LHCb Makefile MinBias.pileup configure cards classes converters display doc examples external modules python readers validation $(DISTDIR)
    437437        @find $(DISTDIR) -depth -name .\* -exec rm -rf {} \;
    438438        @tar -czf $(DISTTAR) $(DISTDIR)

  • external/TrackCovariance/ObsTrk.cc

    r5eb063e ra643a7f  
    1 #include <iostream>
    2 
    31#include <TMath.h>
    42#include <TVector3.h>
    … …  
    75#include <TDecompChol.h>
    86#include <TRandom.h>
    9 
     7#include <iostream>
    108#include "SolGridCov.h"
    119#include "ObsTrk.h"
    12 
    13 using namespace std;
    14 
     10//
     11// Constructors
    1512// x(3) track origin, p(3) track momentum at origin, Q charge, B magnetic field in Tesla
    1613ObsTrk::ObsTrk(TVector3 x, TVector3 p, Double_t Q, Double_t B, SolGridCov *GC)
    1714{
    18   fGC = GC;
    19   fGenX = x;
    20   fGenP = p;
    21   fGenQ = Q;
    22   fB = B;
    23   fGenPar.ResizeTo(5);
    24   fObsPar.ResizeTo(5);
    25   fCov.ResizeTo(5, 5);
    26   fGenPar = XPtoPar(x, p, Q);
    27   fObsPar = GenToObsPar(fGenPar, fGC);
    28   fObsX = ParToX(fObsPar);
    29   fObsP = ParToP(fObsPar);
    30   fObsQ = ParToQ(fObsPar);
    31 }
    32 
     15        fGC = GC;
     16        fGenX = x;
     17        fGenP = p;
     18        fGenQ = Q;
     19        fB = B;
     20        fGenPar.ResizeTo(5);
     21        fGenParACTS.ResizeTo(6);
     22        fGenParILC.ResizeTo(5);
     23        fObsPar.ResizeTo(5);
     24        fObsParACTS.ResizeTo(6);
     25        fObsParILC.ResizeTo(5);
     26        fCov.ResizeTo(5, 5);
     27        fCovACTS.ResizeTo(6, 6);
     28        fCovILC.ResizeTo(5, 5);
     29        fGenPar = XPtoPar(x,p,Q);
     30        fGenParACTS = ParToACTS(fGenPar);
     31        fGenParILC = ParToILC(fGenPar);
     32        /*
     33        std::cout << "ObsTrk::ObsTrk: fGenPar";
     34        for (Int_t i = 0; i < 5; i++)std::cout << fGenPar(i) << ", ";
     35        std::cout << std::endl;
     36        */
     37        fObsPar = GenToObsPar(fGenPar, fGC);
     38        fObsParACTS = ParToACTS(fObsPar);
     39        fObsParILC = ParToILC(fObsPar);
     40        fObsX = ParToX(fObsPar);
     41        fObsP = ParToP(fObsPar);
     42        fObsQ = ParToQ(fObsPar);
     43        fCovACTS = CovToACTS(fCov);
     44        fCovILC = CovToILC(fCov);
     45}
     46//
     47// Destructor
    3348ObsTrk::~ObsTrk()
    3449{
    35 }
    36 
     50        fGenX.Clear();
     51        fGenP.Clear();
     52        fGenPar.Clear();
     53        fGenParACTS.Clear();
     54        fObsX.Clear();
     55        fObsP.Clear();
     56        fObsPar.Clear();
     57        fObsParACTS.Clear();
     58        fCov.Clear();
     59        fCovACTS.Clear();
     60}
    3761TVectorD ObsTrk::XPtoPar(TVector3 x, TVector3 p, Double_t Q)
    3862{
    39   TVectorD Par(5);
    40   // Transverse parameters
    41   Double_t a = -Q * fB * 0.2998; // Units are Tesla, GeV and m
    42   Double_t pt = p.Pt();
    43   Double_t C = a / (2 * pt); // Half curvature
    44 
    45   Double_t r2 = x.Perp2();
    46   Double_t cross = x(0) * p(1) - x(1) * p(0);
    47   Double_t T = TMath::Sqrt(pt * pt - 2 * a * cross + a * a * r2);
    48   Double_t phi0 = TMath::ATan2((p(1) - a * x(0)) / T, (p(0) + a * x(1)) / T); // Phi0
    49   Double_t D; // Impact parameter D
    50   if (pt < 10.0) D = (T - pt) / a;
    51   else D = (-2 * cross + a * r2) / (T + pt);
    52 
    53   Par(0) = D; // Store D
    54   Par(1) = phi0; // Store phi0
    55   Par(2) = C; // Store C
    56   // Longitudinal parameters
    57   Double_t B = C * TMath::Sqrt(TMath::Max(r2 - D * D,0.0) / (1 + 2 * C * D));
    58   Double_t st = TMath::ASin(B) / C;
    59   Double_t ct = p(2) / pt;
    60   Double_t z0 = x(2) - ct * st;
    61 
    62   Par(3) = z0; // Store z0
    63   Par(4) = ct; // Store cot(theta)
    64 
    65   return Par;
    66 }
    67 
     63        //
     64        TVectorD Par(5);
     65        // Transverse parameters
     66        Double_t a = -Q*fB*0.2998;                      // Units are Tesla, GeV and meters
     67        Double_t pt = p.Pt();
     68        Double_t C = a / (2 * pt);                      // Half curvature
     69        //std::cout << "ObsTrk::XPtoPar: fB = " << fB << ", a = " << a << ", pt = " << pt << ", C = " << C << std::endl;
     70        Double_t r2 = x.Perp2();
     71        Double_t cross = x(0)*p(1) - x(1)*p(0);
     72        Double_t T = TMath::Sqrt(pt*pt - 2 * a*cross + a*a*r2);
     73        Double_t phi0 = TMath::ATan2((p(1) - a*x(0)) / T, (p(0) + a*x(1)) / T); // Phi0
     74        Double_t D;                                                     // Impact parameter D
     75        if (pt < 10.0) D = (T - pt) / a;
     76        else D = (-2 * cross + a*r2) / (T + pt);
     77        //
     78        Par(0) = D;             // Store D
     79        Par(1) = phi0;  // Store phi0
     80        Par(2) = C;             // Store C
     81        //Longitudinal parameters
     82        Double_t B = C*TMath::Sqrt(TMath::Max(r2 - D*D,0.0) / (1 + 2 * C*D));
     83        Double_t st = TMath::ASin(B) / C;
     84        Double_t ct = p(2) / pt;
     85        Double_t z0 = x(2) - ct*st;
     86        //
     87        Par(3) = z0;            // Store z0
     88        Par(4) = ct;            // Store cot(theta)
     89        //
     90        return Par;
     91}
     92//
    6893TVector3 ObsTrk::ParToX(TVectorD Par)
    6994{
    70   Double_t D    = Par(0);
    71   Double_t phi0 = Par(1);
    72   Double_t z0   = Par(3);
    73 
    74   TVector3 Xval;
    75   Xval(0) = -D * TMath::Sin(phi0);
    76   Xval(1) =  D * TMath::Cos(phi0);
    77   Xval(2) =  z0;
    78 
    79   return Xval;
    80 }
    81 
     95        Double_t D    = Par(0);
     96        Double_t phi0 = Par(1);
     97        Double_t z0   = Par(3);
     98        //
     99        TVector3 Xval;
     100        Xval(0) = -D*TMath::Sin(phi0);
     101        Xval(1) =  D*TMath::Cos(phi0); 
     102        Xval(2) =  z0;
     103        //
     104        return Xval;
     105}
     106//
    82107TVector3 ObsTrk::ParToP(TVectorD Par)
    83108{
    84   Double_t C    = Par(2);
    85   Double_t phi0 = Par(1);
    86   Double_t ct   = Par(4);
    87   //
    88   TVector3 Pval;
    89   Double_t pt = fB * 0.2998 / TMath::Abs(2 * C);
    90   Pval(0) = pt * TMath::Cos(phi0);
    91   Pval(1) = pt * TMath::Sin(phi0);
    92   Pval(2) = pt * ct;
    93   //
    94   return Pval;
    95 }
     109        Double_t C    = Par(2);
     110        Double_t phi0 = Par(1);
     111        Double_t ct   = Par(4);
     112        //
     113        TVector3 Pval;
     114        Double_t pt = fB*0.2998 / TMath::Abs(2 * C);
     115        Pval(0) = pt*TMath::Cos(phi0);
     116        Pval(1) = pt*TMath::Sin(phi0);
     117        Pval(2) = pt*ct;
     118        //
     119        return Pval;
     120}
     121//
    96122
    97123Double_t ObsTrk::ParToQ(TVectorD Par)
    98124{
    99   return TMath::Sign(1.0, -Par(2));
    100 }
    101 
     125        return TMath::Sign(1.0, -Par(2));
     126}
     127//
    102128TVectorD ObsTrk::GenToObsPar(TVectorD gPar, SolGridCov *GC)
    103129{
    104   TVector3 p = ParToP(gPar);
    105   Double_t pt = p.Pt();
    106   Double_t tanTh = 1.0 / TMath::Abs(gPar(4));
    107   Double_t angd = TMath::ATan(tanTh) * 180. / TMath::Pi();
    108   // Check ranges
    109   Double_t minPt = GC->GetMinPt ();
    110   if (pt < minPt) cout << "Warning ObsTrk::GenToObsPar: pt " << pt << " is below grid range of " << minPt << endl;
    111   Double_t maxPt = GC->GetMaxPt();
    112   if (pt > maxPt) cout << "Warning ObsTrk::GenToObsPar: pt " << pt << " is above grid range of " << maxPt << endl;
    113   Double_t minAn = GC->GetMinAng();
    114   if (angd < minAn) cout << "Warning ObsTrk::GenToObsPar: angle " << angd
    115     << " is below grid range of " << minAn << endl;
    116   Double_t maxAn = GC->GetMaxAng();
    117   if (angd > maxAn) cout << "Warning ObsTrk::GenToObsPar: angle " << angd
    118     << " is above grid range of " << maxAn << endl;
    119   TMatrixDSym Cov = GC->GetCov(pt, angd);
    120   fCov = Cov;
    121   // Now do Choleski decomposition and random number extraction, with appropriate stabilization
    122   TMatrixDSym CvN = Cov;
    123   TMatrixDSym DCv(5); DCv.Zero();
    124   TMatrixDSym DCvInv(5); DCvInv.Zero();
    125   for (Int_t id = 0; id < 5; id++)
    126   {
    127     Double_t dVal = TMath::Sqrt(Cov(id, id));
    128     DCv   (id, id) = dVal;
    129     DCvInv(id, id) = 1.0 / dVal;
    130   }
    131   CvN.Similarity(DCvInv); // Normalize diagonal to 1
    132   TDecompChol Chl(CvN);
    133   Bool_t OK = Chl.Decompose(); // Choleski decomposition of normalized matrix
    134   TMatrixD U = Chl.GetU(); // Get Upper triangular matrix
    135   TMatrixD Ut(TMatrixD::kTransposed, U); // Transposed of U (lower triangular)
    136   TVectorD r(5);
    137   for (Int_t i = 0; i < 5; i++) r(i) = gRandom->Gaus(0.0, 1.0); // Array of normal random numbers
    138   TVectorD oPar = gPar + DCv * (Ut * r); // Observed parameter vector
    139 
    140   return oPar;
    141 }
     130        TVector3 p = ParToP(gPar);
     131        Double_t pt = p.Pt();
     132        Double_t tanTh = 1.0 / TMath::Abs(gPar(4));
     133        Double_t angd = TMath::ATan(tanTh)*180. / TMath::Pi();
     134        //
     135        // Check ranges
     136        Double_t minPt = GC->GetMinPt ();
     137        if (pt < minPt) std::cout << "Warning ObsTrk::GenToObsPar: pt " << pt << " is below grid range of " << minPt << std::endl;
     138        Double_t maxPt = GC->GetMaxPt();
     139        if (pt > maxPt) std::cout << "Warning ObsTrk::GenToObsPar: pt " << pt << " is above grid range of " << maxPt << std::endl;
     140        Double_t minAn = GC->GetMinAng();
     141        if (angd < minAn) std::cout << "Warning ObsTrk::GenToObsPar: angle " << angd
     142                << " is below grid range of " << minAn << std::endl;
     143        Double_t maxAn = GC->GetMaxAng();
     144        if (angd > maxAn) std::cout << "Warning ObsTrk::GenToObsPar: angle " << angd
     145                << " is above grid range of " << maxAn << std::endl;
     146        //
     147        TMatrixDSym Cov = GC->GetCov(pt, angd);
     148        fCov = Cov;
     149        //
     150        // Now do Choleski decomposition and random number extraction, with appropriate stabilization
     151        //
     152        TMatrixDSym CvN = Cov;
     153        TMatrixDSym DCv(5); DCv.Zero();
     154        TMatrixDSym DCvInv(5); DCvInv.Zero();
     155        for (Int_t id = 0; id < 5; id++)
     156        {
     157                Double_t dVal = TMath::Sqrt(Cov(id, id));
     158                DCv   (id, id) = dVal;
     159                DCvInv(id, id) = 1.0 / dVal;
     160        }
     161        CvN.Similarity(DCvInv);                 // Normalize diagonal to 1
     162        TDecompChol Chl(CvN);
     163        Bool_t OK = Chl.Decompose();            // Choleski decomposition of normalized matrix
     164        TMatrixD U = Chl.GetU();                        // Get Upper triangular matrix
     165        TMatrixD Ut(TMatrixD::kTransposed, U); // Transposed of U (lower triangular)
     166        TVectorD r(5);
     167        for (Int_t i = 0; i < 5; i++)r(i) = gRandom->Gaus(0.0, 1.0);            // Array of normal random numbers
     168        TVectorD oPar = gPar + DCv*(Ut*r);      // Observed parameter vector
     169        //
     170        return oPar;
     171}
     172// Parameter conversion to ACTS format
     173TVectorD ObsTrk::ParToACTS(TVectorD Par)
     174{
     175        TVectorD pACTS(6);      // Return vector
     176        //
     177        Double_t b = -0.29988*fB / 2.;
     178        pACTS(0) = 1000*Par(0);         // D from m to mm
     179        pACTS(1) = 1000 * Par(3);       // z0 from m to mm
     180        pACTS(2) = Par(1);                      // Phi0 is unchanged
     181        pACTS(3) = TMath::ATan(1.0 / Par(4)) + TMath::PiOver2();                // Theta in [0, pi] range
     182        pACTS(4) = Par(2) / (b*TMath::Sqrt(1 + Par(4)*Par(4)));         // q/p in GeV
     183        pACTS(5) = 0.0;                         // Time: currently undefined
     184        //
     185        return pACTS;
     186}
     187// Covariance conversion to ACTS format
     188TMatrixDSym ObsTrk::CovToACTS(TMatrixDSym Cov)
     189{
     190        TMatrixDSym cACTS(6); cACTS.Zero();
     191        Double_t b = -0.29988*fB / 2.;
     192        //
     193        // Fill derivative matrix
     194        TMatrixD A(5, 5);       A.Zero();
     195        Double_t ct = fGenPar(4);       // cot(theta)
     196        Double_t C = fGenPar(2);                // half curvature
     197        A(0, 0) = 1000.;                // D-D  conversion to mm
     198        A(1, 2) = 1.0;          // phi0-phi0
     199        A(2, 4) = 1.0/(TMath::Sqrt(1.0 + ct*ct) * b);   // q/p-C
     200        A(3, 1) = 1000.;                // z0-z0 conversion to mm
     201        A(4, 3) = -1.0 / (1.0 + ct*ct); // theta - cot(theta)
     202        A(4, 4) = -C*ct / (b*pow(1.0 + ct*ct,3.0/2.0)); // q/p-cot(theta)
     203        //
     204        TMatrixDSym Cv = Cov;
     205        TMatrixD At(5, 5);
     206        At.Transpose(A);
     207        Cv.Similarity(At);
     208        TMatrixDSub(cACTS, 0, 4, 0, 4) = Cv;
     209        cACTS(5, 5) = 0.1;      // Currently undefined: set to arbitrary value to avoid crashes
     210        //
     211        return cACTS;
     212}
     213
     214// Parameter conversion to ILC format
     215TVectorD ObsTrk::ParToILC(TVectorD Par)
     216{
     217        TVectorD pILC(5);       // Return vector
     218        //
     219        pILC(0) = Par(0)*1.0e3;                 // d0 in mm
     220        pILC(1) = Par(1);                               // phi0 is unchanged
     221        pILC(2) = -2 * Par(2)*1.0e-3;   // w in mm^-1
     222        pILC(3) = Par(3)*1.0e3;                 // z0 in mm
     223        pILC(4) = Par(4);                               // tan(lambda) = cot(theta)
     224        //
     225        return pILC;
     226}
     227// Covariance conversion to ILC format
     228TMatrixDSym ObsTrk::CovToILC(TMatrixDSym Cov)
     229{
     230        TMatrixDSym cILC(5); cILC.Zero();
     231        //
     232        // Fill derivative matrix
     233        TMatrixD A(5, 5);       A.Zero();
     234        //
     235        A(0, 0) = 1.0e3;                // D-d0 in mm
     236        A(1, 1) = 1.0;          // phi0-phi0
     237        A(2, 2) = -2.0e-3;      // w-C
     238        A(3, 3) = 1.0e3;                // z0-z0 conversion to mm
     239        A(4, 4) = 1.0;          // tan(lambda) - cot(theta)
     240        //
     241        TMatrixDSym Cv = Cov;
     242        TMatrixD At(5, 5);
     243        At.Transpose(A);
     244        Cv.Similarity(At);
     245        cILC = Cv;
     246        //
     247        return cILC;
     248}
     249
     250
     251
     252       

  • external/TrackCovariance/ObsTrk.h

    r5eb063e ra643a7f  
     1//
    12#ifndef G__OBSTRK_H
    23#define G__OBSTRK_H
    3 
    44#include <TVector3.h>
    55#include <TVectorD.h>
    66#include <TMatrixDSym.h>
    7 
    8 class SolGridCov;
    9 
     7#include <TDecompChol.h>
     8#include "SolGridCov.h"
     9//
    1010// Class to handle smearing of generated charged particle tracks
    11 
     11//
    1212// Author: F. Bedeschi
    1313//         INFN - Sezione di Pisa, Italy
    1414//
    1515class ObsTrk{
    16   // Class to handle simulation of tracking resolution
    17   // Prefix Obs marks variables after resolution smearing
    18   // Prefix Gen marks variables before resolution smearing
    19 private:
    20   Double_t fB;      // Solenoid magnetic field
    21   SolGridCov *fGC;  // Covariance matrix grid
    22   Double_t fGenQ;   // Generated track charge
    23   Double_t fObsQ;   // Observed  track charge
    24   TVector3 fGenX;   // Generated track origin (x,y,z)
    25   TVector3 fObsX;   // Observed  track origin (x,y,z) @ track min. approach
    26   TVector3 fGenP;   // Generated track momentum at track origin
    27   TVector3 fObsP;   // Observed  track momentum @ track minimum approach
    28   TVectorD fGenPar; // Generated helix track parameters (D, phi0, C, z0, cot(th))
    29   TVectorD fObsPar; // Observed  helix track parameters (D, phi0, C, z0, cot(th))
    30   TMatrixDSym fCov; // INterpolated covariance of track parameters
     16        //
     17        // Class to handle simulation of tracking resolution
     18        // Prefix Obs marks variables after resolution smearing
     19        // Prefix Gen marks variables before resolution smearing
     20        //
     21private:       
     22        Double_t fB;                                            // Solenoid magnetic field
     23        SolGridCov *fGC;                                        // Covariance matrix grid
     24        Double_t fGenQ;                                 // Generated track charge
     25        Double_t fObsQ;                                 // Observed  track charge
     26        TVector3 fGenX;                                 // Generated track origin (x,y,z)
     27        TVector3 fObsX;                                 // Observed  track origin (x,y,z) @ track min. approach
     28        TVector3 fGenP;                                 // Generated track momentum at track origin
     29        TVector3 fObsP;                                 // Observed  track momentum @ track minimum approach
     30        TVectorD fGenPar;                               // Generated helix track parameters (D, phi0, C, z0, cot(th))
     31        TVectorD fGenParACTS;                   // Generated helix track parameters (D, z0, phi0, th, q/p, time
     32        TVectorD fGenParILC;                            // Generated helix track parameters (w, phi0, d0, z0, tan(lambda))
     33        TVectorD fObsPar;                               // Observed  helix track parameters (D, phi0, C, z0, cot(th))
     34        TVectorD fObsParACTS;                   // Observed  helix track parameters (D, z0, phi0, th, q/p, time
     35        TVectorD fObsParILC;                            // Observed  helix track parameters (d0, phi0, w, z0, tan(lambda))
     36        TMatrixDSym fCov;                               // INterpolated covariance of track parameters
     37        TMatrixDSym fCovACTS;                   // Covariance of track parameters in ACTS format
     38                                                                        // (D, z0, phi0, theta, q/p, time)
     39        TMatrixDSym fCovILC;                            // Covariance of track parameters in ILC format
     40                                                                        // (d0, phi0, w, z0, tan(lambda))
     41        //
     42        // Conversion to ACTS parametrization
     43        //
     44        TVectorD ParToACTS(TVectorD Par);               // Parameter conversion
     45        TMatrixDSym CovToACTS(TMatrixDSym Cov); // Covariance
     46        //
     47        // Conversion to ILC parametrization
     48        //
     49        TVectorD ParToILC(TVectorD Par);                // Parameter conversion
     50        TMatrixDSym CovToILC(TMatrixDSym Cov);  // Covariance conversion
     51        //
    3152public:
    32   // x(3) track origin, p(3) track momentum at origin, Q charge, B magnetic field in Tesla
    33   ObsTrk(TVector3 x, TVector3 p, Double_t Q, Double_t B, SolGridCov *GC); // Initialize and generate smeared track
    34   ~ObsTrk();
    35   // Service routines
    36   TVectorD XPtoPar(TVector3 x, TVector3 p, Double_t Q);
    37   TVectorD GenToObsPar(TVectorD gPar, SolGridCov *GC);
    38   TVector3 ParToX(TVectorD Par);
    39   TVector3 ParToP(TVectorD Par);
    40   Double_t ParToQ(TVectorD Par);
    41   // Accessors
    42   // Generator level X, P, Q
    43   Double_t GetGenQ() { return fGenQ; }
    44   TVector3 GetGenX() { return fGenX; }
    45   TVector3 GetGenP() { return fGenP; }
    46   // D, phi0, C, z0, cot(th)
    47   TVectorD GetGenPar() { return fGenPar; }
    48   // Observed level X, P, Q
    49   Double_t GetObsQ() { return fObsQ; }
    50   TVector3 GetObsX() { return fObsX; }
    51   TVector3 GetObsP() { return fObsP; }
    52   // D, phi0, C, z0, cot(th)
    53   TVectorD GetObsPar() { return fObsPar; }
    54   TMatrixDSym GetCov() { return fCov; }
     53        //
     54        // Constructors
     55        // x(3) track origin, p(3) track momentum at origin, Q charge, B magnetic field in Tesla
     56        ObsTrk(TVector3 x, TVector3 p, Double_t Q, Double_t B, SolGridCov *GC); // Initialize and generate smeared track
     57        // Destructor
     58        ~ObsTrk();
     59        //
     60        // Service routines
     61        //
     62        TVectorD XPtoPar(TVector3 x, TVector3 p, Double_t Q);
     63        TVectorD GenToObsPar(TVectorD gPar, SolGridCov *GC);
     64        TVector3 ParToX(TVectorD Par);
     65        TVector3 ParToP(TVectorD Par);
     66        Double_t ParToQ(TVectorD Par);
     67        //
     68        // Accessors
     69        //
     70        // Generator level:
     71        // X, P, Q
     72        Double_t GetGenQ()      { return fGenQ; }
     73        TVector3 GetGenX()      { return fGenX; }
     74        TVector3 GetGenP()      { return fGenP; }
     75        // D, phi0, C, z0, cot(th)
     76        TVectorD GetGenPar()    { return fGenPar; }
     77        // D, z0, phi0, theta, q/p, time
     78        TVectorD GetGenParACTS()        { return fGenParACTS; }
     79        // d0, phi0, w, z0, tan(lambda)
     80        TVectorD GetGenParILC() { return fGenParILC; }
     81        // Observed level X, P, Q
     82        Double_t GetObsQ()      { return fObsQ; }
     83        TVector3 GetObsX()      { return fObsX; }
     84        TVector3 GetObsP()      { return fObsP; }
     85        // D, phi0, C, z0, cot(th)
     86        TVectorD GetObsPar()    { return fObsPar; }
     87        // D, z0, phi0, theta, q/p, time
     88        TVectorD GetObsParACTS()        { return fObsParACTS; }
     89        // d0, phi0, w, z0, tan(lambda)
     90        TVectorD GetObsParILC() { return fObsParILC; }
     91        // Covariances
     92        TMatrixDSym GetCov(){ return fCov; }
     93        TMatrixDSym GetCovACTS(){ return fCovACTS; }
     94        TMatrixDSym GetCovILC(){ return fCovILC; }
    5595};
    5696

  • modules/Efficiency.cc

    r5eb063e ra643a7f  
    7878  fItInputArray = fInputArray->MakeIterator();
    7979
     80  // switch to compute efficiency based on momentum vector eta, phi
     81  fUseMomentumVector = GetBool("UseMomentumVector", false);
     82
    8083  // create output array
    8184
    … …  
    104107    eta = candidatePosition.Eta();
    105108    phi = candidatePosition.Phi();
     109
     110    if (fUseMomentumVector){
     111      eta = candidateMomentum.Eta();
     112      phi = candidateMomentum.Phi();
     113    }
     114
    106115    pt = candidateMomentum.Pt();
    107116    e = candidateMomentum.E();

  • modules/Efficiency.h

    r5eb063e ra643a7f  
    5353  TObjArray *fOutputArray; //!
    5454
     55  Double_t fUseMomentumVector; //!
     56
    5557  ClassDef(Efficiency, 1)
    5658};

  • modules/FastJetFinder.cc

    r5eb063e ra643a7f  
    489489      fastjet::Filter trimmer(fastjet::JetDefinition(fastjet::kt_algorithm, fRTrim), fastjet::SelectorPtFractionMin(fPtFracTrim));
    490490      fastjet::PseudoJet trimmed_jet = trimmer(*itOutputList);
    491 
    492       trimmed_jet = join(trimmed_jet.constituents());
    493 
     491     
    494492      candidate->TrimmedP4[0].SetPtEtaPhiM(trimmed_jet.pt(), trimmed_jet.eta(), trimmed_jet.phi(), trimmed_jet.m());
    495493

  • modules/TrackCovariance.cc

    r5eb063e ra643a7f  
    9696  Candidate *candidate, *mother;
    9797  Double_t mass, p, pt, q, ct;
    98   Double_t dd0, ddz, dphi, dct, dp, dpt;
     98  Double_t dd0, ddz, dphi, dct, dp, dpt, dC;
    9999 
    100100
    … …  
    114114    candidate->Momentum.SetVectM(track.GetObsP(), mass);
    115115    candidate->InitialPosition.SetXYZT(track.GetObsX().X(),track.GetObsX().Y(),track.GetObsX().Z(),candidatePosition.T());
     116   
     117    // save full covariance 5x5 matrix internally (D0, phi, Curvature, dz, ctg(theta))
     118    candidate->TrackCovariance = track.GetCov();
    116119
    117120    pt = candidate->Momentum.Pt();
    … …  
    125128   
    126129    candidate->D0 = track.GetObsPar()[0];
     130    candidate->Phi = track.GetObsPar()[1];
     131    candidate->C = track.GetObsPar()[2];
    127132    candidate->DZ = track.GetObsPar()[3];
     133    candidate->CtgTheta = track.GetObsPar()[4];
    128134    candidate->P  = track.GetObsP().Mag();
    129     candidate->CtgTheta = track.GetObsPar()[4];
    130     candidate->Phi = track.GetObsPar()[1];
    131 
    132135    candidate->PT = pt;
    133136    candidate->Charge = q;
    … …  
    139142    dpt       = 2 * TMath::Sqrt( track.GetCov()(2, 2))*pt*pt / (0.2998*fBz);
    140143    dp        = TMath::Sqrt((1.+ct*ct)*dpt*dpt + 4*pt*pt*ct*ct*dct*dct/(1.+ct*ct)/(1.+ct*ct));
     144    dC        = TMath::Sqrt(track.GetCov()(2, 2));
    141145
    142146    candidate->ErrorD0 = dd0;
    143147    candidate->ErrorDZ = ddz;
    144148    candidate->ErrorP = dp;
     149    candidate->ErrorC = dC;
    145150    candidate->ErrorCtgTheta = dct;
    146151    candidate->ErrorPhi = dphi;
    … …  
    148153    //candidate->TrackResolution = dpt / pt;
    149154    candidate->TrackResolution = dp / p;
     155
    150156
    151157    candidate->AddCandidate(mother);

  • validation/DelphesValidation.cpp

    r5eb063e ra643a7f  
    3232#include "TGraph.h"
    3333#include "TGraphErrors.h"
     34#include "TF1.h"
    3435#include "TH1.h"
    3536#include "TH2.h"

  • validation/validation.sh

    r5eb063e ra643a7f  
    3838outputrootdir=report_${cardbase%.*}/root
    3939cardlabel=${cardbase%.*}
    40 version=$(cat VERSION)
     40version=x.y.z
    4141outpdf=$mainoutputdir/${output%.*}.pdf
    4242cardsdir=validation/cards
Note: See TracChangeset for help on using the changeset viewer.