1 #include <algo/inc/jointprob.h>
4 #include <util/inc/helixrep.h>
5 #include <util/inc/projection.h>
10 #include <util/inc/string.h>
17 _MaxD0Significance = 200;
21 temp.push_back(1.01313412);
22 temp.push_back(0.0246350896);
23 temp.push_back(0.102197811);
24 temp.push_back(0.0411203019);
25 temp.push_back(0.0157710761);
27 _ResolutionParameterRphi = temp;
30 temp1.push_back(1.01629865);
31 temp1.push_back(0.0271386635);
32 temp1.push_back(0.0948112309);
33 temp1.push_back(0.0410759225);
34 temp1.push_back(0.0148685882);
36 _ResolutionParameterZ = temp1;
39 temp2.push_back(1.02015948);
40 temp2.push_back(0.0177643765);
41 temp2.push_back(0.144750029);
42 temp2.push_back(0.0288017225);
43 temp2.push_back(0.0237413906);
45 _ResolutionParameter3D = temp2;
47 _ParameterNames.push_back(
"MaxD0Significance");
48 _ParameterNames.push_back(
"MaxD0andZ0");
49 _ParameterNames.push_back(
"ResolutionParameterRphi");
50 _ParameterNames.push_back(
"ResolutionParameterZ");
52 _ParameterValues.push_back(makeString(_MaxD0Significance));
53 _ParameterValues.push_back(makeString(_MaxD0andZ0));
54 _ParameterValues.push_back(
"Default Values");
55 _ParameterValues.push_back(
"Default Values");
56 _ParameterValues.push_back(
"Default Values");
67 return _ParameterNames;
72 _ParameterValues.clear();
73 _ParameterValues.push_back(makeString(_MaxD0Significance));
74 _ParameterValues.push_back(makeString(_MaxD0andZ0));
75 _ParameterValues.push_back(makeString(_ResolutionParameterRphi[0]));
76 _ParameterValues.push_back(makeString(_ResolutionParameterRphi[1]));
77 _ParameterValues.push_back(makeString(_ResolutionParameterRphi[2]));
78 _ParameterValues.push_back(makeString(_ResolutionParameterRphi[3]));
79 _ParameterValues.push_back(makeString(_ResolutionParameterRphi[4]));
80 _ParameterValues.push_back(makeString(_ResolutionParameterZ[0]));
81 _ParameterValues.push_back(makeString(_ResolutionParameterZ[1]));
82 _ParameterValues.push_back(makeString(_ResolutionParameterZ[2]));
83 _ParameterValues.push_back(makeString(_ResolutionParameterZ[3]));
84 _ParameterValues.push_back(makeString(_ResolutionParameterZ[4]));
86 return _ParameterValues;
91 this->badParameter(Parameter);
96 if (Parameter ==
"MaxD0Significance")
98 _MaxD0Significance = Value;
101 if (Parameter ==
"MaxD0andZ0")
106 this->badParameter(Parameter);
111 if (Parameter ==
"ResolutionParameterRphi")
113 _ResolutionParameterRphi = *(std::vector<double>*) Value;
116 if (Parameter ==
"ResolutionParameterZ")
118 _ResolutionParameterZ = *(std::vector<double>*) Value;
121 if (Parameter ==
"ResolutionParameter3D")
123 _ResolutionParameter3D = *(std::vector<double>*) Value;
126 else this->badParameter(Parameter);
132 std::map<Projection,double> ResultMap;
135 double totalprod[3] = {1,1,1};
136 double jprob[3] = {0,0,0};
137 int ntraks[3] = {0,0,0};
140 double maxdz = _MaxD0andZ0;
141 double maxd0sig = _MaxD0Significance;
142 double significancecompare = 0;
148 for (std::vector<Track*>::const_iterator iTrack= (MyJet->
tracks().begin()); iTrack != (MyJet->
tracks().end()) ;++iTrack)
152 if( fabs( (*iTrack)->helixRep().d0() )< maxdz && fabs( (*iTrack)->helixRep().z0() ) < maxdz )
156 for(j = 0; j<3; j++ )
162 significancecompare = fabs((*iTrack)->significance(RPhi));
167 significancecompare = fabs((*iTrack)->significance( Z ));
173 significancecompare = fabs((*iTrack)->significance(ThreeD));
178 if( significancecompare < maxd0sig )
182 totalprod[j] *= probparam( significancecompare ,j , maxd0sig )/ probparam( 0 , j, maxd0sig );
193 for( j = 0; j<3; j++ )
201 for( k =0; k < ntraks[j]; k++ )
204 if( k > 0 ) fact *= k;
209 sigma += pow( ( -log( totalprod[j] ) ), k ) / fact;
213 jprob[j] = totalprod[j] * sigma;
219 ResultMap[RPhi] = jprob[0];
220 ResultMap[Z] = jprob[1];
221 ResultMap[ThreeD] = jprob[2];
226 double JointProb::probparam(
double parameter,
int thecoord,
double maxd0sig )
const
230 double resolutionparameters[3][5];
231 if (_ResolutionParameterRphi.size() != 5 || _ResolutionParameterZ.size() != 5 || _ResolutionParameter3D.size() != 5)
232 std::cerr <<
"Warning jointprob.cpp:229 Parameters of wrong length" << std::endl;
233 for(
int iii=0; iii<5; iii++)
235 resolutionparameters[0][iii] = _ResolutionParameterRphi[iii];
236 resolutionparameters[1][iii] = _ResolutionParameterZ[iii];
237 resolutionparameters[2][iii] = _ResolutionParameter3D[iii];
245 prob = erfc( parameter / (sqrt(
double(2) ) * resolutionparameters[thecoord][0] ))
246 -erfc( maxd0sig / ( sqrt(
double(2) ) * resolutionparameters[thecoord][0] ) );
250 prob += resolutionparameters[thecoord][1]* ( exp(- resolutionparameters[thecoord][2] * parameter )
251 -exp(- resolutionparameters[thecoord][2] * maxd0sig ) )
252 +resolutionparameters[thecoord][3]* ( exp(- resolutionparameters[thecoord][4] * parameter )
253 -exp(- resolutionparameters[thecoord][4] * maxd0sig ) );
262 prob = exp(- ( parameter * parameter ) / ( resolutionparameters[thecoord][0] * resolutionparameters[thecoord][0] *
double ( 2 ) ) )
263 -exp(- ( maxd0sig * maxd0sig ) / ( resolutionparameters[thecoord][0] * resolutionparameters[thecoord][0] *
double ( 2 ) ) );
265 prob += resolutionparameters[thecoord][1]* ( ( 1 + resolutionparameters[thecoord][2] * parameter )
266 *exp ( - resolutionparameters[thecoord][2] * parameter )
267 - ( 1 + resolutionparameters[thecoord][2] * maxd0sig )
268 *exp ( - resolutionparameters[thecoord][2] * maxd0sig ))
269 +resolutionparameters[thecoord][3]* ( ( 1 + resolutionparameters[thecoord][4] * parameter )
270 *exp ( - resolutionparameters[thecoord][4] * parameter )
271 - ( 1 + resolutionparameters[thecoord][4] * maxd0sig )
272 *exp ( - resolutionparameters[thecoord][4] * maxd0sig ));
std::map< Projection, double > calculateFor(Jet *MyJet) const
Run the algorithm on a jet.
double prob(double ChiSquared, double DegreesOfFreedom)
Calculate probability from Chi Squared and Degrees of freedom.
void setStringParameter(const string &Parameter, const string &Value)
Set String Parameter.
const std::vector< Track * > & tracks() const
Tracks.
void setPointerParameter(const string &Parameter, void *Value)
Set Pointer Parameter.
void setDoubleParameter(const string &Parameter, const double Value)
Set Double Parameter.
std::vector< string > parameterValues() const
Parameter Values.
string name() const
Name.
std::vector< string > parameterNames() const
Parameter Names.