001 #include "GaudiKernel/ToolFactory.h"
002 #include "GaudiKernel/SmartDataPtr.h"
003 #include "GaudiKernel/IDataProviderSvc.h"
004
005
006 #include "TrkParameters/MeasuredPerigee.h"
007
008
009 #include "STACOTools/StacoCombineToolNG01.h"
010
011 StacoCombineToolNG01::StacoCombineToolNG01(const std::string& t,
012 const std::string& n,
013 const IInterface* p ):AthAlgTool(t,n,p)
014 {
015 declareInterface<IStacoCombineTool>(this);
016
017 declareProperty("PreSelOpt1" , m_PreSelOpt1 = 1 ) ;
018 declareProperty("PreSelOpt1DeltaR" , m_PreSelOpt1DeltaR = 0.5 ) ;
019
020 declareProperty("PreSelOpt2" , m_PreSelOpt2 = 0 ) ;
021 declareProperty("PreSelOpt2DeltaZ0CTE" , m_PreSelOpt2DeltaZ0CTE = 800. ) ;
022 declareProperty("PreSelOpt2DeltaZ0SLOPE" , m_PreSelOpt2DeltaZ0SLOPE = 0. ) ;
023 declareProperty("PreSelOpt2DeltaEtaCTE" , m_PreSelOpt2DeltaEtaCTE = 0.25 ) ;
024 declareProperty("PreSelOpt2DeltaEtaSLOPE" , m_PreSelOpt2DeltaEtaSLOPE = 0. ) ;
025
026 }
027
028 StacoCombineToolNG01::~StacoCombineToolNG01(){}
029
030
031 StatusCode StacoCombineToolNG01::initialize(){
032
033 msg(MSG::INFO) << "Initialisation started " << endreq;
034
035 StatusCode sc = AthAlgTool::initialize();
036 if ( sc.isFailure() ) {
037 msg(MSG::FATAL) << " AthAlgTool::initialize() failed" << endreq;
038 return( StatusCode::FAILURE );
039 }
040
041 msg(MSG::INFO) << "================================" << endreq;
042 msg(MSG::INFO) << "=Proprieties are " << endreq;
043 msg(MSG::INFO) << "= PreSelOpt1 " << m_PreSelOpt1 << endreq;
044 msg(MSG::INFO) << "= PreSelOpt1DeltaR " << m_PreSelOpt1DeltaR << endreq;
045 msg(MSG::INFO) << "= " << endreq;
046 msg(MSG::INFO) << "= PreSelOpt2 " << m_PreSelOpt2 << endreq;
047 msg(MSG::INFO) << "= PreSelOpt2DeltaZ0CTE " << m_PreSelOpt2DeltaZ0CTE << endreq;
048 msg(MSG::INFO) << "= PreSelOpt2DeltaZ0SLOPE " << m_PreSelOpt2DeltaZ0SLOPE << endreq;
049 msg(MSG::INFO) << "= PreSelOpt2DeltaEtaCTE " << m_PreSelOpt2DeltaEtaCTE << endreq;
050 msg(MSG::INFO) << "= PreSelOpt2DeltaEtaSLOPE " << m_PreSelOpt2DeltaEtaSLOPE << endreq;
051 msg(MSG::INFO) << "================================" << endreq;
052
053
054 msg(MSG::INFO) << "Initialisation ended " << endreq;
055 return StatusCode::SUCCESS;
056
057 }
058
059 StatusCode StacoCombineToolNG01::finalize(){return StatusCode::SUCCESS;}
060
061
062
063 bool StacoCombineToolNG01::TheCombIdMu(
064 const Trk::MeasuredPerigee* pMeasuredPerigeeID,
065 const Trk::MeasuredPerigee* pMeasuredPerigeeMS,
066 double& Xi2,
067 double& ThePTInv
068 ){
069
070 ThePTInv = 1e20;
071 double Xi2Mu = 0. ;
072 double Xi2ID = 0. ;
073 double ParCB[5],covCB[15];
074
075 if (!TheCombIdMu(pMeasuredPerigeeID, pMeasuredPerigeeMS,
076 Xi2, Xi2Mu, Xi2ID, ParCB, covCB
077 )) return false;
078
079 ThePTInv =ParCB[4]/sin(ParCB[3]);
080
081 return true;
082
083 }
084
085 bool StacoCombineToolNG01::TheCombIdMu(
086 const Trk::MeasuredPerigee* pMeasuredPerigeeID,
087 const Trk::MeasuredPerigee* pMeasuredPerigeeMS,
088 double& Xi2,double& Xi2ID,double& Xi2Mu,
089 double* ParCB,double* covCB
090 ){
091
092 Xi2 = 1e20 ;
093 Xi2Mu = 0. ;
094 Xi2ID = 0. ;
095
096 if ( pMeasuredPerigeeID == 0 ) return false;
097
098 if ( pMeasuredPerigeeMS == 0 ) return false;
099
100 if ( !Preselection(pMeasuredPerigeeID,pMeasuredPerigeeMS) ) return false;
101
102 const Trk::ErrorMatrix ErrorMatrix_ID = pMeasuredPerigeeID->localErrorMatrix();
103 const Trk::WeightMatrix* pWeightMatrix_ID = ErrorMatrix_ID.weightPointer();
104 if ( pWeightMatrix_ID == 0 ){
105 msg(MSG::INFO) << " pWeightMatrix_ID null " << endreq;
106 return false ;
107 }
108 const Trk::CovarianceMatrix* pCovarianceMatrix_ID = ErrorMatrix_ID.covariancePointer();
109 if ( pCovarianceMatrix_ID == 0 ){
110 msg(MSG::INFO) << " pCovarianceMatrix_ID null " << endreq;
111 return false ;
112 }
113
114 const Trk::ErrorMatrix ErrorMatrix_MS = pMeasuredPerigeeMS->localErrorMatrix();
115 const Trk::WeightMatrix* pWeightMatrix_MS = ErrorMatrix_MS.weightPointer();
116 if ( pWeightMatrix_MS == 0 ){
117 msg(MSG::INFO) << " pWeightMatrix_MS null " << endreq;
118 return false ;
119 }
120 const Trk::CovarianceMatrix* pCovarianceMatrix_MS = ErrorMatrix_MS.covariancePointer();
121 if ( pCovarianceMatrix_MS == 0 ){
122 msg(MSG::INFO) << " pCovarianceMatrix_MS null " << endreq;
123 return false ;
124 }
125
126 int ifail;
127
128 Trk::WeightMatrix WeightMatrix_CB = (*pWeightMatrix_ID) + (*pWeightMatrix_MS) ;
129 Trk::CovarianceMatrix CovarianceMatrix_CB(WeightMatrix_CB.inverse(ifail)) ;
130 if (ifail){
131 msg(MSG::INFO) << " Inversion of WeightMatrix_CB failed " << endreq;
132 return false ;
133 }
134
135 Trk::CovarianceMatrix SummOfCovarianceMatrix = (*pCovarianceMatrix_ID) + (*pCovarianceMatrix_MS) ;
136 Trk::WeightMatrix InvSummOfCovarianceMatrix(SummOfCovarianceMatrix.inverse(ifail)) ;
137 if (ifail){
138 msg(MSG::INFO) << " Inversion of SummOfCovarianceMatrix failed " << endreq;
139 return false ;
140 }
141
142 HepVector Parameters_CB = CovarianceMatrix_CB * ( (*pWeightMatrix_ID) * pMeasuredPerigeeID->parameters() + (*pWeightMatrix_MS) * pMeasuredPerigeeMS->parameters() ) ;
143
144 HepVector DiffParameters(5) ;
145
146 DiffParameters[Trk::d0] = pMeasuredPerigeeID->parameters()[Trk::d0] - pMeasuredPerigeeMS->parameters()[Trk::d0] ;
147 DiffParameters[Trk::z0] = pMeasuredPerigeeID->parameters()[Trk::z0] - pMeasuredPerigeeMS->parameters()[Trk::z0] ;
148 DiffParameters[Trk::phi] = pMeasuredPerigeeID->parameters()[Trk::phi] - pMeasuredPerigeeMS->parameters()[Trk::phi] ;
149 DiffParameters[Trk::theta] = pMeasuredPerigeeID->parameters()[Trk::theta] - pMeasuredPerigeeMS->parameters()[Trk::theta] ;
150 DiffParameters[Trk::qOverP] = pMeasuredPerigeeID->parameters()[Trk::qOverP] - pMeasuredPerigeeMS->parameters()[Trk::qOverP] ;
151 if(DiffParameters[Trk::phi]>M_PI) DiffParameters[Trk::phi] -= 2.*M_PI; else if(DiffParameters[Trk::phi]<-M_PI) DiffParameters[Trk::phi] +=2.*M_PI;
152 Xi2 = InvSummOfCovarianceMatrix.similarity(DiffParameters) ;
153 Xi2 = Xi2/5. ;
154
155 DiffParameters[Trk::d0] = Parameters_CB[Trk::d0] - pMeasuredPerigeeID->parameters()[Trk::d0] ;
156 DiffParameters[Trk::z0] = Parameters_CB[Trk::z0] - pMeasuredPerigeeID->parameters()[Trk::z0] ;
157 DiffParameters[Trk::phi] = Parameters_CB[Trk::phi] - pMeasuredPerigeeID->parameters()[Trk::phi] ;
158 DiffParameters[Trk::theta] = Parameters_CB[Trk::theta] - pMeasuredPerigeeID->parameters()[Trk::theta] ;
159 DiffParameters[Trk::qOverP] = Parameters_CB[Trk::qOverP] - pMeasuredPerigeeID->parameters()[Trk::qOverP] ;
160 if(DiffParameters[Trk::phi]>M_PI) DiffParameters[Trk::phi] -= 2.*M_PI; else if(DiffParameters[Trk::phi]<-M_PI) DiffParameters[Trk::phi] +=2.*M_PI;
161 Xi2ID = (*pWeightMatrix_ID).similarity(DiffParameters) ;
162 Xi2ID = Xi2ID/5. ;
163
164 DiffParameters[Trk::d0] = Parameters_CB[Trk::d0] - pMeasuredPerigeeMS->parameters()[Trk::d0] ;
165 DiffParameters[Trk::z0] = Parameters_CB[Trk::z0] - pMeasuredPerigeeMS->parameters()[Trk::z0] ;
166 DiffParameters[Trk::phi] = Parameters_CB[Trk::phi] - pMeasuredPerigeeMS->parameters()[Trk::phi] ;
167 DiffParameters[Trk::theta] = Parameters_CB[Trk::theta] - pMeasuredPerigeeMS->parameters()[Trk::theta] ;
168 DiffParameters[Trk::qOverP] = Parameters_CB[Trk::qOverP] - pMeasuredPerigeeMS->parameters()[Trk::qOverP] ;
169 if(DiffParameters[Trk::phi]>M_PI) DiffParameters[Trk::phi] -= 2.*M_PI; else if(DiffParameters[Trk::phi]<-M_PI) DiffParameters[Trk::phi] +=2.*M_PI;
170 Xi2Mu = (*pWeightMatrix_MS).similarity(DiffParameters) ;
171 Xi2ID = Xi2ID/5. ;
172
173 ParCB[0] = Parameters_CB[Trk::d0];
174 ParCB[1] = Parameters_CB[Trk::z0];
175 ParCB[2] = Parameters_CB[Trk::phi];
176 ParCB[3] = Parameters_CB[Trk::theta];
177 ParCB[4] = Parameters_CB[Trk::qOverP];
178
179 covCB[0] = CovarianceMatrix_CB[Trk::d0][Trk::d0];
180 covCB[1] = CovarianceMatrix_CB[Trk::z0][Trk::d0];
181 covCB[2] = CovarianceMatrix_CB[Trk::z0][Trk::z0];
182 covCB[3] = CovarianceMatrix_CB[Trk::phi][Trk::d0];
183 covCB[4] = CovarianceMatrix_CB[Trk::phi][Trk::z0];
184 covCB[5] = CovarianceMatrix_CB[Trk::phi][Trk::phi];
185 covCB[6] = CovarianceMatrix_CB[Trk::theta][Trk::d0];
186 covCB[7] = CovarianceMatrix_CB[Trk::theta][Trk::z0];
187 covCB[8] = CovarianceMatrix_CB[Trk::theta][Trk::phi];
188 covCB[9] = CovarianceMatrix_CB[Trk::theta][Trk::theta];
189 covCB[10] = CovarianceMatrix_CB[Trk::qOverP][Trk::d0];
190 covCB[11] = CovarianceMatrix_CB[Trk::qOverP][Trk::z0];
191 covCB[12] = CovarianceMatrix_CB[Trk::qOverP][Trk::phi];
192 covCB[13] = CovarianceMatrix_CB[Trk::qOverP][Trk::theta];
193 covCB[14] = CovarianceMatrix_CB[Trk::qOverP][Trk::qOverP];
194
195 return true;
196
197 }
198
199
200
201 void StacoCombineToolNG01::GetPar_TrackParticle(
202 const Trk::MeasuredPerigee* pMeasuredPerigee,
203 double* Param,double* CovMat
204 ){
205
206 for(int i=0; i<5; i++) {Param[i]=0.;}
207 for(int i=0; i<15; i++) {CovMat[i]=0.;}
208
209 if ( pMeasuredPerigee == 0 ) return ;
210
211
212
213
214 Param[0] = pMeasuredPerigee->parameters()[Trk::d0];
215 Param[1] = pMeasuredPerigee->parameters()[Trk::z0];
216 Param[2] = pMeasuredPerigee->parameters()[Trk::phi];
217 Param[3] = pMeasuredPerigee->parameters()[Trk::theta];
218 Param[4] = pMeasuredPerigee->parameters()[Trk::qOverP];
219
220
221 const Trk::ErrorMatrix& ermatrix=pMeasuredPerigee->localErrorMatrix();
222
223 CovMat[0] = ermatrix.covValue(Trk::d0,Trk::d0);
224 CovMat[1] = ermatrix.covValue(Trk::z0,Trk::d0);
225 CovMat[2] = ermatrix.covValue(Trk::z0,Trk::z0);
226 CovMat[3] = ermatrix.covValue(Trk::phi,Trk::d0);
227 CovMat[4] = ermatrix.covValue(Trk::phi,Trk::z0);
228 CovMat[5] = ermatrix.covValue(Trk::phi,Trk::phi);
229 CovMat[6] = ermatrix.covValue(Trk::theta,Trk::d0);
230 CovMat[7] = ermatrix.covValue(Trk::theta,Trk::z0);
231 CovMat[8] = ermatrix.covValue(Trk::theta,Trk::phi);
232 CovMat[9] = ermatrix.covValue(Trk::theta,Trk::theta);
233 CovMat[10] = ermatrix.covValue(Trk::qOverP,Trk::d0);
234 CovMat[11] = ermatrix.covValue(Trk::qOverP,Trk::z0);
235 CovMat[12] = ermatrix.covValue(Trk::qOverP,Trk::phi);
236 CovMat[13] = ermatrix.covValue(Trk::qOverP,Trk::theta);
237 CovMat[14] = ermatrix.covValue(Trk::qOverP,Trk::qOverP);
238
239 }
240
241
242
243
244 bool StacoCombineToolNG01::TheCombIdMu(double* ,double* ,double* , double& ,double& ,double& ){return false ;}
245 bool StacoCombineToolNG01::Invert55 (double* ,double* ){return false ;}
246 void StacoCombineToolNG01::HelToPar (double* ,double* ,double* ){}
247 void StacoCombineToolNG01::GetPar_TrackParticle(const Trk::MeasuredPerigee* ,double* ){}
248
249
250
251
252
253 bool StacoCombineToolNG01::Preselection(
254 const Trk::MeasuredPerigee* pMeasuredPerigeeID,
255 const Trk::MeasuredPerigee* pMeasuredPerigeeMS
256 ){
257
258 if ( m_PreSelOpt1 == 1 ){
259 if ( !Preselection01(pMeasuredPerigeeID,pMeasuredPerigeeMS) ) return false;
260 }
261 if ( m_PreSelOpt2 == 1 ){
262 if ( !Preselection02(pMeasuredPerigeeID,pMeasuredPerigeeMS) ) return false;
263 }
264
265 return true ;
266 }
267
268 bool StacoCombineToolNG01::Preselection01(
269 const Trk::MeasuredPerigee* pMeasuredPerigeeID,
270 const Trk::MeasuredPerigee* pMeasuredPerigeeMS
271 ){
272
273 double EtaID = asinh(1./tan(pMeasuredPerigeeID->parameters()[Trk::theta]));
274 double EtaMS = asinh(1./tan(pMeasuredPerigeeMS->parameters()[Trk::theta]));
275 double dEta = fabs(EtaID-EtaMS);
276
277 double PhiID = pMeasuredPerigeeID->parameters()[Trk::phi] ;
278 double PhiMS = pMeasuredPerigeeMS->parameters()[Trk::phi] ;
279 double dPhi = fabs(PhiID-PhiMS) ;
280 if (dPhi > M_PI) dPhi = 2.*M_PI - dPhi;
281
282 double DeltaR2 = dEta*dEta+dPhi*dPhi;
283
284 double DeltaRCut = m_PreSelOpt1DeltaR ;
285 if (DeltaR2>DeltaRCut*DeltaRCut) return false;
286
287 return true ;
288 }
289
290 bool StacoCombineToolNG01::Preselection02(
291 const Trk::MeasuredPerigee* pMeasuredPerigeeID,
292 const Trk::MeasuredPerigee* pMeasuredPerigeeMS
293 ){
294
295 double TheInvPID = fabs(1./pMeasuredPerigeeID->parameters()[Trk::qOverP]) ;
296
297 double Z0ID = pMeasuredPerigeeID->parameters()[Trk::z0] ;
298 double Z0MS = pMeasuredPerigeeMS->parameters()[Trk::z0] ;
299 double dZ0 = fabs( Z0ID - Z0MS );
300
301 double DeltaZ0Cut = m_PreSelOpt2DeltaZ0CTE + TheInvPID * m_PreSelOpt2DeltaZ0SLOPE ;
302 if (dZ0>DeltaZ0Cut) return false;
303
304 double EtaID = asinh(1./tan(pMeasuredPerigeeID->parameters()[Trk::theta]));
305 double EtaMS = asinh(1./tan(pMeasuredPerigeeMS->parameters()[Trk::theta]));
306 double dEta = fabs( EtaID - EtaMS );
307
308 double DeltaEtaCut = m_PreSelOpt2DeltaEtaCTE + TheInvPID * m_PreSelOpt2DeltaEtaSLOPE ;
309 if (dEta>DeltaEtaCut) return false;
310
311 return true ;
312 }
| Due to the LXR bug, the updates fail sometimes to remove references to deleted files. The Saturday's full rebuilds fix these problems |
|
This page was automatically generated by the
LXR engine.
|
|