00001
00002 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00003
00004 #include "Alignment/CommonAlignment/interface/Alignable.h"
00005 #include "Alignment/CommonAlignment/interface/AlignableDetOrUnitPtr.h"
00006 #include "Alignment/CommonAlignmentParametrization/interface/CompositeAlignmentDerivativesExtractor.h"
00007 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
00008
00009 #include "Alignment/CommonAlignmentParametrization/interface/CompositeAlignmentParameters.h"
00010
00011
00012
00013 CompositeAlignmentParameters::
00014 CompositeAlignmentParameters(const AlgebraicVector& par, const AlgebraicSymMatrix& cov,
00015 const Components& comp) :
00016 AlignmentParameters(0,par,cov) ,
00017 theComponents(comp)
00018 {}
00019
00020
00021
00022 CompositeAlignmentParameters::
00023 CompositeAlignmentParameters(const AlgebraicVector& par, const AlgebraicSymMatrix& cov,
00024 const Components& comp, const AlignableDetToAlignableMap& alimap,
00025 const Aliposmap& aliposmap, const Alilenmap& alilenmap) :
00026 AlignmentParameters(0,par,cov) ,
00027 theComponents(comp) ,
00028 theAlignableDetToAlignableMap(alimap),
00029 theAliposmap(aliposmap),
00030 theAlilenmap(alilenmap)
00031 {}
00032
00033
00034
00035 CompositeAlignmentParameters::
00036 CompositeAlignmentParameters(const DataContainer& data,
00037 const Components& comp, const AlignableDetToAlignableMap& alimap,
00038 const Aliposmap& aliposmap, const Alilenmap& alilenmap) :
00039 AlignmentParameters(0,data) ,
00040 theComponents(comp) ,
00041 theAlignableDetToAlignableMap(alimap),
00042 theAliposmap(aliposmap),
00043 theAlilenmap(alilenmap)
00044 {}
00045
00046
00047
00048 CompositeAlignmentParameters::~CompositeAlignmentParameters()
00049 {}
00050
00051
00052
00053 CompositeAlignmentParameters*
00054 CompositeAlignmentParameters::clone( const AlgebraicVector& par,
00055 const AlgebraicSymMatrix& cov) const
00056 {
00057 CompositeAlignmentParameters* cap =
00058 new CompositeAlignmentParameters(par,cov,components());
00059
00060 if ( userVariables() )
00061 cap->setUserVariables(userVariables()->clone());
00062
00063 return cap;
00064 }
00065
00066
00067
00068 CompositeAlignmentParameters*
00069 CompositeAlignmentParameters::cloneFromSelected( const AlgebraicVector& par,
00070 const AlgebraicSymMatrix& cov) const
00071 {
00072 return clone(par,cov);
00073 }
00074
00075
00076
00077 CompositeAlignmentParameters*
00078 CompositeAlignmentParameters::clone( const AlgebraicVector& par,
00079 const AlgebraicSymMatrix& cov,
00080 const AlignableDetToAlignableMap& alimap,
00081 const Aliposmap& aliposmap,
00082 const Alilenmap& alilenmap ) const
00083 {
00084 CompositeAlignmentParameters* cap =
00085 new CompositeAlignmentParameters(par,cov,components(),alimap,aliposmap,alilenmap);
00086
00087 if ( userVariables() )
00088 cap->setUserVariables(userVariables()->clone());
00089
00090 return cap;
00091 }
00092
00093
00094
00095 CompositeAlignmentParameters*
00096 CompositeAlignmentParameters::cloneFromSelected( const AlgebraicVector& par,
00097 const AlgebraicSymMatrix& cov,
00098 const AlignableDetToAlignableMap& alimap,
00099 const Aliposmap& aliposmap,
00100 const Alilenmap& alilenmap) const
00101 {
00102 return clone(par,cov,alimap,aliposmap,alilenmap);
00103 }
00104
00105
00106
00107 CompositeAlignmentParameters::Components
00108 CompositeAlignmentParameters::components() const
00109 {
00110 return theComponents;
00111 }
00112
00113
00114
00115
00116 AlgebraicMatrix
00117 CompositeAlignmentParameters::derivatives( const std::vector<TrajectoryStateOnSurface>& tsosvec,
00118 const std::vector<AlignableDet*>& alidetvec ) const
00119 {
00120 std::vector<AlignableDetOrUnitPtr> detOrUnits;
00121 this->convert(alidetvec, detOrUnits);
00122
00123 return this->derivatives(tsosvec, detOrUnits);
00124 }
00125
00126 AlgebraicMatrix
00127 CompositeAlignmentParameters::derivatives( const std::vector<TrajectoryStateOnSurface>& tsosvec,
00128 const std::vector<AlignableDetOrUnitPtr>& alidetvec ) const
00129 {
00130 std::vector<Alignable*> alivec;
00131 for (std::vector<AlignableDetOrUnitPtr>::const_iterator it=alidetvec.begin(); it!=alidetvec.end(); ++it)
00132 alivec.push_back(alignableFromAlignableDet(*it));
00133
00134 CompositeAlignmentDerivativesExtractor extractor(alivec,alidetvec,tsosvec);
00135 return extractor.derivatives();
00136 }
00137
00138
00139 AlgebraicVector
00140 CompositeAlignmentParameters::correctionTerm( const std::vector<TrajectoryStateOnSurface>& tsosvec,
00141 const std::vector<AlignableDet*>& alidetvec) const
00142 {
00143 std::vector<AlignableDetOrUnitPtr> detOrUnits;
00144 this->convert(alidetvec, detOrUnits);
00145
00146 return this->correctionTerm(tsosvec, detOrUnits);
00147 }
00148
00149
00150 AlgebraicVector
00151 CompositeAlignmentParameters::correctionTerm( const std::vector<TrajectoryStateOnSurface>& tsosvec,
00152 const std::vector<AlignableDetOrUnitPtr>& alidetvec) const
00153 {
00154 std::vector<Alignable*> alivec;
00155 for (std::vector<AlignableDetOrUnitPtr>::const_iterator it=alidetvec.begin(); it!=alidetvec.end(); ++it )
00156 alivec.push_back(alignableFromAlignableDet(*it));
00157
00158 CompositeAlignmentDerivativesExtractor extractor(alivec,alidetvec,tsosvec);
00159 return extractor.correctionTerm();
00160 }
00161
00162
00163
00164 AlgebraicMatrix
00165 CompositeAlignmentParameters::selectedDerivatives( const std::vector<TrajectoryStateOnSurface>& tsosvec,
00166 const std::vector<AlignableDet*>& alidetvec) const
00167 {
00168 return derivatives(tsosvec,alidetvec);
00169 }
00170
00171
00172 AlgebraicMatrix
00173 CompositeAlignmentParameters::selectedDerivatives( const std::vector<TrajectoryStateOnSurface>& tsosvec,
00174 const std::vector<AlignableDetOrUnitPtr>& alidetvec) const
00175 {
00176 return derivatives(tsosvec,alidetvec);
00177 }
00178
00179
00180
00181 AlgebraicMatrix
00182 CompositeAlignmentParameters::derivatives( const TrajectoryStateOnSurface &tsos,
00183 const AlignableDetOrUnitPtr &alidet) const
00184 {
00185 std::vector<TrajectoryStateOnSurface> tsosvec;
00186 std::vector<AlignableDetOrUnitPtr> alidetvec;
00187 tsosvec.push_back(tsos);
00188 alidetvec.push_back(alidet);
00189 return derivatives(tsosvec,alidetvec);
00190 }
00191
00192
00193
00194 AlgebraicMatrix
00195 CompositeAlignmentParameters::selectedDerivatives( const TrajectoryStateOnSurface &tsos,
00196 const AlignableDetOrUnitPtr &alidet ) const
00197 {
00198 return derivatives(tsos,alidet);
00199 }
00200
00201
00202
00203
00204
00205 AlgebraicMatrix
00206 CompositeAlignmentParameters::derivativesLegacy( const std::vector<TrajectoryStateOnSurface> &tsosvec,
00207 const std::vector<AlignableDet*>& alidetvec ) const
00208 {
00209
00210 if (alidetvec.size() != tsosvec.size())
00211 {
00212 edm::LogError("BadArgument") << " Inconsistent length of argument vectors! ";
00213 AlgebraicMatrix selderiv(1,0);
00214 return selderiv;
00215 }
00216
00217 std::vector<AlgebraicMatrix> vecderiv;
00218 int nparam=0;
00219
00220 std::vector<TrajectoryStateOnSurface>::const_iterator itsos=tsosvec.begin();
00221 for( std::vector<AlignableDet*>::const_iterator it=alidetvec.begin();
00222 it!=alidetvec.end(); ++it, ++itsos )
00223 {
00224 AlignableDet* ad = (*it);
00225 Alignable* ali = alignableFromAlignableDet(ad);
00226 AlignmentParameters* ap = ali->alignmentParameters();
00227 AlgebraicMatrix thisselderiv = ap->selectedDerivatives(*itsos,ad);
00228 vecderiv.push_back(thisselderiv);
00229 nparam += thisselderiv.num_row();
00230 }
00231
00232 int ipos=1;
00233 AlgebraicMatrix selderiv(nparam,2);
00234 for ( std::vector<AlgebraicMatrix>::const_iterator imat=vecderiv.begin();
00235 imat!=vecderiv.end(); ++imat )
00236 {
00237 AlgebraicMatrix thisselderiv=(*imat);
00238 int npar=thisselderiv.num_row();
00239 selderiv.sub(ipos,1,thisselderiv);
00240 ipos += npar;
00241 }
00242
00243 return selderiv;
00244 }
00245
00246
00247
00248
00249 AlgebraicMatrix
00250 CompositeAlignmentParameters::selectedDerivativesLegacy( const std::vector<TrajectoryStateOnSurface> &tsosvec,
00251 const std::vector<AlignableDet*>& alidetvec ) const
00252 {
00253 return derivativesLegacy(tsosvec,alidetvec);
00254 }
00255
00256
00257
00258
00259 AlgebraicMatrix
00260 CompositeAlignmentParameters::derivativesLegacy( const TrajectoryStateOnSurface& tsos,
00261 AlignableDet* alidet ) const
00262 {
00263 std::vector<TrajectoryStateOnSurface> tsosvec;
00264 std::vector<AlignableDet*> alidetvec;
00265 tsosvec.push_back(tsos);
00266 alidetvec.push_back(alidet);
00267 return derivativesLegacy(tsosvec,alidetvec);
00268
00269 }
00270
00271
00272
00273
00274 AlgebraicMatrix
00275 CompositeAlignmentParameters::selectedDerivativesLegacy( const TrajectoryStateOnSurface& tsos,
00276 AlignableDet* alidet ) const
00277 {
00278 return derivativesLegacy(tsos,alidet);
00279 }
00280
00281
00282
00283
00284 Alignable*
00285 CompositeAlignmentParameters::alignableFromAlignableDet(AlignableDetOrUnitPtr adet) const
00286 {
00287
00288 AlignableDetToAlignableMap::const_iterator iali =
00289 theAlignableDetToAlignableMap.find(adet);
00290 if ( iali!=theAlignableDetToAlignableMap.end() ) return (*iali).second;
00291 else return 0;
00292
00293 }
00294
00295
00296
00297 AlgebraicVector
00298 CompositeAlignmentParameters::parameterSubset( const std::vector<Alignable*>& vec ) const
00299 {
00300 const std::vector< Alignable* > sel = extractAlignables( vec );
00301
00302 const unsigned int nali = sel.size();
00303 int ndim = 0;
00304
00305 std::vector<int> posvec;
00306 std::vector<int> lenvec;
00307
00308 posvec.reserve( nali );
00309 lenvec.reserve( nali );
00310
00311
00312 if ( !extractPositionAndLength( sel, posvec, lenvec, ndim ) ) return AlgebraicVector();
00313
00314
00315 AlgebraicVector result( ndim );
00316
00317 int resi = 0;
00318 for ( unsigned int iali = 0; iali < nali; ++iali )
00319 {
00320 int posi = posvec[iali];
00321 int leni = lenvec[iali];
00322
00323 for ( int ir = 0; ir < leni; ++ir )
00324 result[resi+ir] = theData->parameters()[posi-1+ir];
00325
00326 resi += leni;
00327 }
00328
00329 return result;
00330 }
00331
00332
00333
00334
00335 AlgebraicSymMatrix
00336 CompositeAlignmentParameters::covarianceSubset( const std::vector<Alignable*>& vec ) const
00337 {
00338 const std::vector< Alignable* > sel = extractAlignables( vec );
00339
00340 const unsigned int nali = sel.size();
00341 int ndim = 0;
00342
00343 std::vector<int> posvec;
00344 std::vector<int> lenvec;
00345
00346 posvec.reserve( nali );
00347 lenvec.reserve( nali );
00348
00349
00350
00351 if ( !extractPositionAndLength( sel, posvec, lenvec, ndim ) ) return AlgebraicSymMatrix();
00352
00353
00354 AlgebraicSymMatrix result( ndim );
00355
00356 int resi = 0;
00357 for ( unsigned int iali = 0; iali < nali; ++iali )
00358 {
00359 int posi = posvec[iali];
00360 int leni = lenvec[iali];
00361
00362 int resj = 0;
00363 for ( unsigned int jali = 0; jali <= iali; ++jali )
00364 {
00365 int posj = posvec[jali];
00366 int lenj = lenvec[jali];
00367
00368 for ( int ir = 0; ir < leni; ++ir )
00369 for ( int ic = 0; ic < lenj; ++ic )
00370 result[resi+ir][resj+ic] = theData->covariance()[posi-1+ir][posj-1+ic];
00371
00372 resj += lenj;
00373 }
00374 resi += leni;
00375 }
00376
00377 return result;
00378 }
00379
00380
00381
00382
00383 AlgebraicMatrix
00384 CompositeAlignmentParameters::covarianceSubset( const std::vector<Alignable*>& veci,
00385 const std::vector<Alignable*>& vecj ) const
00386 {
00387 const std::vector< Alignable* > seli = extractAlignables( veci );
00388 const std::vector< Alignable* > selj = extractAlignables( vecj );
00389
00390 int ndimi=0;
00391 int ndimj=0;
00392
00393 std::vector<int> posveci;
00394 std::vector<int> lenveci;
00395 std::vector<int> posvecj;
00396 std::vector<int> lenvecj;
00397
00398 posveci.reserve( seli.size() );
00399 lenveci.reserve( seli.size() );
00400 posvecj.reserve( selj.size() );
00401 lenvecj.reserve( selj.size() );
00402
00403
00404
00405 if ( !extractPositionAndLength( seli, posveci, lenveci, ndimi ) ) return AlgebraicSymMatrix();
00406
00407 if ( !extractPositionAndLength( selj, posvecj, lenvecj, ndimj ) ) return AlgebraicSymMatrix();
00408
00409
00410 AlgebraicMatrix result( ndimi, ndimj );
00411
00412 int resi = 0;
00413 for ( unsigned int iali = 0; iali < seli.size(); ++iali )
00414 {
00415 int posi = posveci[iali];
00416 int leni = lenveci[iali];
00417
00418 int resj = 0;
00419 for ( unsigned int jali = 0; jali < selj.size(); ++jali )
00420 {
00421 int posj = posvecj[jali];
00422 int lenj = lenvecj[jali];
00423
00424 for ( int ir = 0; ir < leni; ++ir )
00425 for ( int ic = 0; ic < lenj; ++ic )
00426 result[resi+ir][resj+ic] = theData->covariance()[posi-1+ir][posj-1+ic];
00427
00428 resj += lenj;
00429 }
00430 resi += leni;
00431 }
00432
00433 return result;
00434 }
00435
00436
00437
00438
00439 bool
00440 CompositeAlignmentParameters::extractPositionAndLength( const std::vector<Alignable*>& alignables,
00441 std::vector<int>& posvec,
00442 std::vector<int>& lenvec,
00443 int& length ) const
00444 {
00445 length = 0;
00446
00447 for ( std::vector<Alignable*>::const_iterator it = alignables.begin(); it != alignables.end(); ++it )
00448 {
00449
00450 if ( std::find( theComponents.begin(), theComponents.end(), *it ) == theComponents.end() )
00451 {
00452 edm::LogError( "NotFound" ) << "@SUB=CompositeAlignmentParameters::extractPositionAndLength"
00453 << "Alignable not found in components!";
00454 return false;
00455 }
00456
00457
00458 Aliposmap::const_iterator iposmap = theAliposmap.find( *it );
00459 Alilenmap::const_iterator ilenmap = theAlilenmap.find( *it );
00460 if ( iposmap == theAliposmap.end() || ilenmap == theAlilenmap.end() )
00461 {
00462 edm::LogError( "NotFound" ) << "@SUB=CompositeAlignmentParameters::extractPositionAndLength"
00463 << "position/length not found for Alignable in maps!";
00464 return false;
00465 }
00466 posvec.push_back( iposmap->second );
00467 lenvec.push_back( ilenmap->second );
00468 length += ilenmap->second;
00469 }
00470
00471 return true;
00472 }
00473
00474
00475
00476 std::vector< Alignable* >
00477 CompositeAlignmentParameters::extractAlignables( const std::vector< Alignable* >& alignables ) const
00478 {
00479 std::vector< Alignable* > result;
00480
00481 std::vector< Alignable* >::const_iterator itA, itEnd;
00482 for ( itA = alignables.begin(), itEnd = alignables.end(); itA != itEnd; ++itA )
00483 {
00484 if ( std::find( result.begin(), result.end(), *itA ) == result.end() ) result.push_back( *itA );
00485 }
00486
00487 return result;
00488 }
00489
00490
00491
00492 void CompositeAlignmentParameters::convert(const std::vector<AlignableDet*> &input,
00493 std::vector<AlignableDetOrUnitPtr> &output) const
00494 {
00495 output.clear();
00496 output.reserve(input.size());
00497
00498 std::vector<AlignableDet*>::const_iterator it, itEnd;
00499 for (it = input.begin(), itEnd = input.end(); it != itEnd; ++it)
00500 output.push_back(AlignableDetOrUnitPtr(*it));
00501 }