MillePedeAlignmentAlgorithm/plugins/MillePedeAlignmentAlgorithm.h

0001 #ifndef Alignment_MillePedeAlignmentAlgorithm_MillePedeAlignmentAlgorithm_h
0002 #define Alignment_MillePedeAlignmentAlgorithm_MillePedeAlignmentAlgorithm_h
0003
0004 /// \class MillePedeAlignmentAlgorithm
0005 ///
0006 ///  CMSSW interface to pede: produces pede's binary input and steering file(s)
0007 ///
0008 ///  \author    : Gero Flucke
0009 ///  date       : October 2006
0010 ///  $Revision: 1.36 $
0011 ///  $Date: 2012/08/10 09:01:11 $
0012 ///  (last update by $Author: flucke $)
0013
0014 #include "Alignment/CommonAlignmentAlgorithm/interface/AlignmentAlgorithmPluginFactory.h"
0015 #include "Alignment/CommonAlignmentAlgorithm/interface/AlignmentAlgorithmBase.h"
0016
0017 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0018 #include "FWCore/Framework/interface/ConsumesCollector.h"
0019
0020 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
0021 #include "TrackingTools/TransientTrackingRecHit/interface/TransientTrackingRecHit.h"
0022
0023 #include "Alignment/ReferenceTrajectories/interface/ReferenceTrajectoryBase.h"
0024
0025 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0026 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0027
0028 #include "Geometry/TrackerGeometryBuilder/interface/PixelTopologyMap.h"
0029 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0030 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0031
0032 #include "CondFormats/PCLConfig/interface/AlignPCLThresholdsHG.h"
0033 #include "CondFormats/DataRecord/interface/AlignPCLThresholdsHGRcd.h"
0034
0035 #include <vector>
0036 #include <string>
0037 #include <memory>
0038
0039 class Alignable;
0040 class AlignableTracker;
0041 class AlignableMuon;
0042 class AlignableExtras;
0043
0044 class AlignmentParameters;
0045 class AlignableNavigator;
0046 class AlignableDetOrUnitPtr;
0047 class AlignmentUserVariables;
0048
0049 class AlignmentParameterStore;
0050
0051 class IntegratedCalibrationBase;
0052
0053 class MillePedeMonitor;
0054 class PedeSteerer;
0055 class PedeLabelerBase;
0056 class Mille;
0057 class TrajectoryFactoryBase;
0058
0059 // already from base class - and forward declaration does not work since typedef!
0060 /* class TkFittedLasBeamCollection; */
0061 /* class TsosVectorCollection; */
0062
0063 class MillePedeAlignmentAlgorithm : public AlignmentAlgorithmBase {
0064 public:
0065   /// Constructor
0066   MillePedeAlignmentAlgorithm(const edm::ParameterSet &cfg, edm::ConsumesCollector &iC);
0067
0068   /// Destructor
0069   ~MillePedeAlignmentAlgorithm() override;
0070
0071   /// Called at beginning of job
0072   void initialize(const edm::EventSetup &setup,
0073                   AlignableTracker *tracker,
0074                   AlignableMuon *muon,
0075                   AlignableExtras *extras,
0076                   AlignmentParameterStore *store) override;
0077
0078   /// Returns whether MP supports calibrations
0079   bool supportsCalibrations() override;
0080   /// Pass integrated calibrations to Millepede (they are not owned by Millepede!)
0081   bool addCalibrations(const std::vector<IntegratedCalibrationBase *> &iCals) override;
0082
0083   virtual bool storeThresholds(const int &nRecords,
0084                                const AlignPCLThresholdsHG::threshold_map &thresholdMap,
0085                                const AlignPCLThresholdsHG::param_map &floatMap);
0086
0087   /// Called at end of job
0088   void terminate(const edm::EventSetup &iSetup) override;
0089   /// Called at end of job
0090   void terminate() override;
0091
0092   /// Returns whether MP should process events in the current configuration
0093   bool processesEvents() override;
0094
0095   /// Returns whether MP produced results to be stored
0096   bool storeAlignments() override;
0097
0098   /// Run the algorithm on trajectories and tracks
0099   void run(const edm::EventSetup &setup, const EventInfo &eventInfo) override;
0100
0101   /// called at begin of run
0102   void beginRun(const edm::Run &run, const edm::EventSetup &setup, bool changed) override;
0103
0104   // TODO: This method does NOT match endRun() in base class! Nobody is
0105   //       calling this?
0106   /// Run on run products, e.g. TkLAS
0107   virtual void endRun(const EventInfo &, const EndRunInfo &,
0108                       const edm::EventSetup &);  //override;
0109
0110   // This one will be called since it matches the interface of the base class
0111   void endRun(const EndRunInfo &runInfo, const edm::EventSetup &setup) override;
0112
0113   /// called at begin of luminosity block (resets Mille binary in mille mode)
0114   void beginLuminosityBlock(const edm::EventSetup &) override;
0115
0116   /// called at end of luminosity block
0117   void endLuminosityBlock(const edm::EventSetup &) override;
0118
0119   /// Called in order to pass parameters to alignables for a specific run
0120   /// range in case the algorithm supports run range dependent alignment.
0121   bool setParametersForRunRange(const RunRange &runrange) override;
0122
0123 private:
0124   enum MeasurementDirection { kLocalX = 0, kLocalY };
0125
0126   /// fill mille for a trajectory, returning number of x/y hits ([0,0] if 'bad' trajectory)
0127   std::pair<unsigned int, unsigned int> addReferenceTrajectory(
0128       const edm::EventSetup &setup,
0129       const EventInfo &eventInfo,
0130       const ReferenceTrajectoryBase::ReferenceTrajectoryPtr &refTrajPtr);
0131
0132   /// If hit is usable: callMille for x and (probably) y direction.
0133   /// If globalDerivatives fine: returns 2 if 2D-hit, 1 if 1D-hit, 0 if no Alignable for hit.
0134   /// Returns -1 if any problem (for params cf. globalDerivativesHierarchy)
0135   int addMeasurementData(const edm::EventSetup &setup,
0136                          const EventInfo &eventInfo,
0137                          const ReferenceTrajectoryBase::ReferenceTrajectoryPtr &refTrajPtr,
0138                          unsigned int iHit,
0139                          AlignmentParameters *&params);
0140
0141   /// Add global data (labels, derivatives) to GBL trajectory
0142   /// Returns -1 if any problem (for params cf. globalDerivativesHierarchy)
0143   int addGlobalData(const edm::EventSetup &setup,
0144                     const EventInfo &eventInfo,
0145                     const ReferenceTrajectoryBase::ReferenceTrajectoryPtr &refTrajPtr,
0146                     unsigned int iHit,
0147                     gbl::GblPoint &gblPoint);
0148
0149   /// Increase hit counting of MillePedeVariables behind each parVec[i]
0150   /// (and also for parameters higher in hierarchy),
0151   /// assuming 'parVec' and 'validHitVecY' to be parallel.
0152   /// Returns number of valid y-hits.
0153   unsigned int addHitCount(const std::vector<AlignmentParameters *> &parVec,
0154                            const std::vector<bool> &validHitVecY) const;
0155
0156   /// adds data from reference trajectory from a specific Hit
0157   template <typename CovarianceMatrix, typename ResidualMatrix, typename LocalDerivativeMatrix>
0158   void addRefTrackData2D(const ReferenceTrajectoryBase::ReferenceTrajectoryPtr &refTrajPtr,
0159                          unsigned int iTrajHit,
0160                          Eigen::MatrixBase<CovarianceMatrix> &aHitCovarianceM,
0161                          Eigen::MatrixBase<ResidualMatrix> &aHitResidualsM,
0162                          Eigen::MatrixBase<LocalDerivativeMatrix> &aLocalDerivativesM);
0163
0164   /// adds data for virtual measurements from reference trajectory
0165   void addVirtualMeas(const ReferenceTrajectoryBase::ReferenceTrajectoryPtr &refTrajPtr, unsigned int iVirtualMeas);
0166
0167   /// adds data for a specific virtual measurement from reference trajectory
0168   template <typename CovarianceMatrix, typename ResidualMatrix, typename LocalDerivativeMatrix>
0169   void addRefTrackVirtualMeas1D(const ReferenceTrajectoryBase::ReferenceTrajectoryPtr &refTrajPtr,
0170                                 unsigned int iVirtualMeas,
0171                                 Eigen::MatrixBase<CovarianceMatrix> &aHitCovarianceM,
0172                                 Eigen::MatrixBase<ResidualMatrix> &aHitResidualsM,
0173                                 Eigen::MatrixBase<LocalDerivativeMatrix> &aLocalDerivativesM);
0174
0175   /// recursively adding derivatives and labels, false if problems
0176   bool globalDerivativesHierarchy(const EventInfo &eventInfo,
0177                                   const TrajectoryStateOnSurface &tsos,
0178                                   Alignable *ali,
0179                                   const AlignableDetOrUnitPtr &alidet,
0180                                   std::vector<float> &globalDerivativesX,
0181                                   std::vector<float> &globalDerivativesY,
0182                                   std::vector<int> &globalLabels,
0183                                   AlignmentParameters *&lowestParams) const;
0184
0185   /// recursively adding derivatives (double) and labels, false if problems
0186   bool globalDerivativesHierarchy(const EventInfo &eventInfo,
0187                                   const TrajectoryStateOnSurface &tsos,
0188                                   Alignable *ali,
0189                                   const AlignableDetOrUnitPtr &alidet,
0190                                   std::vector<double> &globalDerivativesX,
0191                                   std::vector<double> &globalDerivativesY,
0192                                   std::vector<int> &globalLabels,
0193                                   AlignmentParameters *&lowestParams) const;
0194
0195   /// adding derivatives from integrated calibrations
0196   void globalDerivativesCalibration(const TransientTrackingRecHit::ConstRecHitPointer &recHit,
0197                                     const TrajectoryStateOnSurface &tsos,
0198                                     const edm::EventSetup &setup,
0199                                     const EventInfo &eventInfo,
0200                                     std::vector<float> &globalDerivativesX,
0201                                     std::vector<float> &globalDerivativesY,
0202                                     std::vector<int> &globalLabels) const;
0203
0204   /// calls callMille1D or callMille2D
0205   int callMille(const ReferenceTrajectoryBase::ReferenceTrajectoryPtr &refTrajPtr,
0206                 unsigned int iTrajHit,
0207                 const std::vector<int> &globalLabels,
0208                 const std::vector<float> &globalDerivativesX,
0209                 const std::vector<float> &globalDerivativesY);
0210
0211   /// calls Mille for 1D hits
0212   int callMille1D(const ReferenceTrajectoryBase::ReferenceTrajectoryPtr &refTrajPtr,
0213                   unsigned int iTrajHit,
0214                   const std::vector<int> &globalLabels,
0215                   const std::vector<float> &globalDerivativesX);
0216
0217   /// calls Mille for x and possibly y component of hit,
0218   /// y is skipped for non-real 2D (e.g. SiStripRecHit2D),
0219   /// for TID/TEC first diagonalises if correlation is larger than configurable
0220   int callMille2D(const ReferenceTrajectoryBase::ReferenceTrajectoryPtr &refTrajPtr,
0221                   unsigned int iTrajHit,
0222                   const std::vector<int> &globalLabels,
0223                   const std::vector<float> &globalDerivativesx,
0224                   const std::vector<float> &globalDerivativesy);
0225
0226   template <typename CovarianceMatrix,
0227             typename LocalDerivativeMatrix,
0228             typename ResidualMatrix,
0229             typename GlobalDerivativeMatrix>
0230   void diagonalize(Eigen::MatrixBase<CovarianceMatrix> &aHitCovarianceM,
0231                    Eigen::MatrixBase<LocalDerivativeMatrix> &aLocalDerivativesM,
0232                    Eigen::MatrixBase<ResidualMatrix> &aHitResidualsM,
0233                    Eigen::MatrixBase<GlobalDerivativeMatrix> &aGlobalDerivativesM) const;
0234
0235   // deals with the non matrix format of theFloatBufferX ...
0236   template <typename GlobalDerivativeMatrix>
0237   void makeGlobDerivMatrix(const std::vector<float> &globalDerivativesx,
0238                            const std::vector<float> &globalDerivativesy,
0239                            Eigen::MatrixBase<GlobalDerivativeMatrix> &aGlobalDerivativesM);
0240
0241   //   void callMille(const ReferenceTrajectoryBase::ReferenceTrajectoryPtr &refTrajPtr,
0242   //               unsigned int iTrajHit, MeasurementDirection xOrY,
0243   //               const std::vector<float> &globalDerivatives, const std::vector<int> &globalLabels);
0244   /// true if hit belongs to 2D detector (currently tracker specific)
0245   bool is2D(const TransientTrackingRecHit::ConstRecHitPointer &recHit) const;
0246
0247   /// read pede input defined by 'psetName', flag to create/not create MillePedeVariables
0248   bool readFromPede(const edm::ParameterSet &mprespset, bool setUserVars, const RunRange &runrange);
0249   bool areEmptyParams(const align::Alignables &alignables) const;
0250   unsigned int doIO(int loop) const;
0251   /// add MillePedeVariables for each AlignmentParameters (exception if no parameters...)
0252   void buildUserVariables(const align::Alignables &alignables) const;
0253
0254   /// Generates list of files to read, given the list and dir from the configuration.
0255   /// This will automatically expand formatting directives, if they appear.
0256   std::vector<std::string> getExistingFormattedFiles(const std::vector<std::string> &plainFiles,
0257                                                      const std::string &theDir);
0258
0259   void addLaserData(const EventInfo &eventInfo,
0260                     const TkFittedLasBeamCollection &tkLasBeams,
0261                     const TsosVectorCollection &tkLasBeamTsoses);
0262   void addLasBeam(const EventInfo &eventInfo,
0263                   const TkFittedLasBeam &lasBeam,
0264                   const std::vector<TrajectoryStateOnSurface> &tsoses);
0265
0266   /// add measurement data from PXB survey
0267   void addPxbSurvey(const edm::ParameterSet &pxbSurveyCfg);
0268
0269   //
0270   bool areIOVsSpecified() const;
0271
0272   //--------------------------------------------------------
0273   // Data members
0274   //--------------------------------------------------------
0275
0276   const edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> topoToken_;
0277   const edm::ESGetToken<AlignPCLThresholdsHG, AlignPCLThresholdsHGRcd> aliThrToken_;
0278   const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> geomToken_;
0279
0280   enum EModeBit { myMilleBit = 1 << 0, myPedeRunBit = 1 << 1, myPedeSteerBit = 1 << 2, myPedeReadBit = 1 << 3 };
0281   unsigned int decodeMode(const std::string &mode) const;
0282   bool isMode(unsigned int testMode) const { return (theMode & testMode); }
0283   bool addHitStatistics(int fromLoop, const std::string &outFile, const std::vector<std::string> &inFiles) const;
0284   bool addHits(const align::Alignables &alis, const std::vector<AlignmentUserVariables *> &mpVars) const;
0285
0286   edm::ParameterSet theConfig;
0287   unsigned int theMode;
0288   std::string theDir;  /// directory for all kind of files
0289   AlignmentParameterStore *theAlignmentParameterStore;
0290   align::Alignables theAlignables;
0291   std::unique_ptr<AlignableNavigator> theAlignableNavigator;
0292   std::unique_ptr<MillePedeMonitor> theMonitor;
0293   std::unique_ptr<Mille> theMille;
0294   std::shared_ptr<PedeLabelerBase> thePedeLabels;
0295   std::unique_ptr<PedeSteerer> thePedeSteer;
0296   std::unique_ptr<TrajectoryFactoryBase> theTrajectoryFactory;
0297   std::vector<IntegratedCalibrationBase *> theCalibrations;
0298   std::shared_ptr<AlignPCLThresholdsHG> theThresholds;
0299   std::shared_ptr<PixelTopologyMap> pixelTopologyMap;
0300   unsigned int theMinNumHits;
0301   double theMaximalCor2D;  /// maximal correlation allowed for 2D hit in TID/TEC.
0302                            /// If larger, the 2D measurement gets diagonalized!!!
0303   const align::RunNumber firstIOV_;
0304   const bool ignoreFirstIOVCheck_;
0305   const bool enableAlignableUpdates_;
0306   int theLastWrittenIov;  // keeping track for output trees...
0307   std::vector<float> theFloatBufferX;
0308   std::vector<float> theFloatBufferY;
0309   std::vector<int> theIntBuffer;
0310   bool theDoSurveyPixelBarrel;
0311   // CHK for GBL
0312   std::unique_ptr<gbl::MilleBinary> theBinary;
0313   bool theGblDoubleBinary;
0314
0315   const bool runAtPCL_;
0316   const bool ignoreHitsWithoutGlobalDerivatives_;
0317   const bool skipGlobalPositionRcdCheck_;
0318
0319   const align::RunRanges uniqueRunRanges_;
0320   const bool enforceSingleIOVInput_;
0321   std::vector<align::RunNumber> cachedRuns_;
0322   align::RunNumber lastProcessedRun_;
0323 };
0324
0325 DEFINE_EDM_PLUGIN(AlignmentAlgorithmPluginFactory, MillePedeAlignmentAlgorithm, "MillePedeAlignmentAlgorithm");
0326 #endif