CMS 3D CMS Logo
List of all members | Public Member Functions | Private Member Functions | Private Attributes
GlobalTrackerMuonAlignment Class Reference

#include <Alignment/GlobalTrackerMuonAlignment/src/GlobalTrackerMuonAlignment.cc>

Inheritance diagram for GlobalTrackerMuonAlignment:
edm::one::EDAnalyzer<> edm::one::EDAnalyzerBase edm::EDConsumerBase

Public Member Functions

void analyzeTrackTrack (const edm::Event &, const edm::EventSetup &)
 
void analyzeTrackTrajectory (const edm::Event &, const edm::EventSetup &)
 
void bookHist ()
 
double CLHEP_dot (const CLHEP::HepVector &a, const CLHEP::HepVector &b)
 
void debugTrackHit (const std::string, reco::TrackRef)
 
void debugTrackHit (const std::string, reco::TransientTrack &)
 
void debugTrajectory (const std::string, Trajectory &)
 
void debugTrajectorySOS (const std::string, TrajectoryStateOnSurface &)
 
void debugTrajectorySOSv (const std::string, TrajectoryStateOnSurface)
 
void fitHist ()
 
 GlobalTrackerMuonAlignment (const edm::ParameterSet &)
 
void gradientGlobal (GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, AlgebraicSymMatrix66 &)
 
void gradientGlobalAlg (GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, AlgebraicSymMatrix66 &)
 
void gradientLocal (GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, CLHEP::HepSymMatrix &, CLHEP::HepMatrix &, CLHEP::HepVector &, AlgebraicVector4 &)
 
void misalignMuon (GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &)
 
void misalignMuonL (GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, AlgebraicVector4 &, TrajectoryStateOnSurface &, TrajectoryStateOnSurface &)
 
void muonFitter (reco::TrackRef, reco::TransientTrack &, PropagationDirection, TrajectoryStateOnSurface &)
 
void trackFitter (reco::TrackRef, reco::TransientTrack &, PropagationDirection, TrajectoryStateOnSurface &)
 
void writeGlPosRcd (CLHEP::HepVector &d3)
 
 ~GlobalTrackerMuonAlignment () override
 
- Public Member Functions inherited from edm::one::EDAnalyzer<>
 EDAnalyzer ()=default
 
 EDAnalyzer (const EDAnalyzer &)=delete
 
SerialTaskQueue * globalLuminosityBlocksQueue () final
 
SerialTaskQueue * globalRunsQueue () final
 
const EDAnalyzer & operator= (const EDAnalyzer &)=delete
 
bool wantsGlobalLuminosityBlocks () const final
 
bool wantsGlobalRuns () const final
 
bool wantsInputProcessBlocks () const final
 
bool wantsProcessBlocks () const final
 
- Public Member Functions inherited from edm::one::EDAnalyzerBase
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
 EDAnalyzerBase ()
 
ModuleDescription const & moduleDescription () const
 
bool wantsStreamLuminosityBlocks () const
 
bool wantsStreamRuns () const
 
 ~EDAnalyzerBase () override
 
- Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfo > consumesInfo () const
 
void convertCurrentProcessAlias (std::string const &processName)
 Convert "@currentProcess" in InputTag process names to the actual current process name. More...
 
 EDConsumerBase ()
 
 EDConsumerBase (EDConsumerBase const &)=delete
 
 EDConsumerBase (EDConsumerBase &&)=default
 
ESProxyIndex const * esGetTokenIndices (edm::Transition iTrans) const
 
std::vector< ESProxyIndex > const & esGetTokenIndicesVector (edm::Transition iTrans) const
 
std::vector< ESRecordIndex > const & esGetTokenRecordIndicesVector (edm::Transition iTrans) const
 
ProductResolverIndexAndSkipBit indexFrom (EDGetToken, BranchType, TypeID const &) const
 
void itemsMayGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
void itemsToGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
std::vector< ProductResolverIndexAndSkipBit > const & itemsToGetFrom (BranchType iType) const
 
void labelsForToken (EDGetToken iToken, Labels &oLabels) const
 
void modulesWhoseProductsAreConsumed (std::array< std::vector< ModuleDescription const *> *, NumBranchTypes > &modulesAll, std::vector< ModuleProcessName > &modulesInPreviousProcesses, ProductRegistry const &preg, std::map< std::string, ModuleDescription const *> const &labelsToDesc, std::string const &processName) const
 
EDConsumerBase const & operator= (EDConsumerBase const &)=delete
 
EDConsumerBase & operator= (EDConsumerBase &&)=default
 
bool registeredToConsume (ProductResolverIndex, bool, BranchType) const
 
bool registeredToConsumeMany (TypeID const &, BranchType) const
 
void selectInputProcessBlocks (ProductRegistry const &productRegistry, ProcessBlockHelperBase const &processBlockHelperBase)
 
ProductResolverIndexAndSkipBit uncheckedIndexFrom (EDGetToken) const
 
void updateLookup (BranchType iBranchType, ProductResolverIndexHelper const &, bool iPrefetchMayGet)
 
void updateLookup (eventsetup::ESRecordsToProxyIndices const &)
 
virtual ~EDConsumerBase () noexcept(false)
 

Private Member Functions

void analyze (const edm::Event &, const edm::EventSetup &) override
 
void beginJob () override
 
void endJob () override
 

Private Attributes

bool collectionCosmic
 
bool collectionIsolated
 
bool cosmicMuonMode_
 
bool debug_
 
bool defineFitter
 
TFile * file
 
AlgebraicVector3 Gfr
 
const Alignments * globalPositionRcd_
 
edm::InputTag gmuonTags_
 
CLHEP::HepVector Grad
 
CLHEP::HepVector GradL
 
CLHEP::HepMatrix Hess
 
CLHEP::HepMatrix HessL
 
TH1F * histo
 
TH2F * histo101
 
TH2F * histo102
 
TH1F * histo11
 
TH1F * histo12
 
TH1F * histo13
 
TH1F * histo14
 
TH1F * histo15
 
TH1F * histo16
 
TH1F * histo17
 
TH1F * histo18
 
TH1F * histo19
 
TH1F * histo2
 
TH1F * histo20
 
TH1F * histo21
 
TH1F * histo22
 
TH1F * histo23
 
TH1F * histo24
 
TH1F * histo25
 
TH1F * histo26
 
TH1F * histo27
 
TH1F * histo28
 
TH1F * histo29
 
TH1F * histo3
 
TH1F * histo30
 
TH1F * histo31
 
TH1F * histo32
 
TH1F * histo33
 
TH1F * histo34
 
TH1F * histo35
 
TH1F * histo4
 
TH1F * histo5
 
TH1F * histo6
 
TH1F * histo7
 
TH1F * histo8
 
AlgebraicSymMatrix33 Inf
 
bool isolatedMuonMode_
 
std::vector< AlignTransform >::const_iterator iteratorEcalRcd
 
std::vector< AlignTransform >::const_iterator iteratorHcalRcd
 
std::vector< AlignTransform >::const_iterator iteratorMuonRcd
 
std::vector< AlignTransform >::const_iterator iteratorTrackerRcd
 
edm::ESGetToken< Alignments, GlobalPositionRcd > m_globalPosToken
 
edm::ESGetToken< MagneticField, IdealMagneticFieldRecord > m_MagFieldToken
 
edm::ESGetToken< Propagator, TrackingComponentsRecord > m_propToken
 
edm::ESGetToken< GlobalTrackingGeometry, GlobalTrackingGeometryRecord > m_TkGeometryToken
 
edm::ESGetToken< TransientTrackingRecHitBuilder, TransientRecHitRecord > m_ttrhBuilderToken
 
const MagneticField * magneticField_
 
CLHEP::HepVector MuGlAngle
 
CLHEP::HepVector MuGlShift
 
edm::InputTag muonTags_
 
MuonTransientTrackingRecHitBuilder * MuRHBuilder
 
std::string MuSelect
 
int N_event
 
int N_track
 
ofstream OutGlobalTxt
 
string propagator_
 
bool refitMuon_
 
bool refitTrack_
 
string rootOutFile_
 
edm::InputTag smuonTags_
 
KFTrajectoryFitter * theFitter
 
KFTrajectoryFitter * theFitterOp
 
KFTrajectorySmoother * theSmoother
 
KFTrajectorySmoother * theSmootherOp
 
edm::ESHandle< GlobalTrackingGeometry > theTrackingGeometry
 
const GlobalTrackingGeometry * trackingGeometry_
 
edm::InputTag trackTags_
 
const TransientTrackingRecHitBuilder * TTRHBuilder
 
string txtOutFile_
 
edm::ESWatcher< GlobalPositionRcd > watchGlobalPositionRcd_
 
edm::ESWatcher< IdealMagneticFieldRecord > watchMagneticFieldRecord_
 
edm::ESWatcher< TrackingComponentsRecord > watchTrackingComponentsRecord_
 
edm::ESWatcher< GlobalTrackingGeometryRecord > watchTrackingGeometry_
 
bool writeDB_
 

Additional Inherited Members

- Public Types inherited from edm::one::EDAnalyzerBase
typedef EDAnalyzerBase ModuleType
 
- Public Types inherited from edm::EDConsumerBase
typedef ProductLabels Labels
 
- Static Public Member Functions inherited from edm::one::EDAnalyzerBase
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &descriptions)
 
- Protected Member Functions inherited from edm::EDConsumerBase
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > consumes (edm::InputTag const &tag)
 
template<BranchType B = InEvent>
EDConsumerBaseAdaptor< B > consumes (edm::InputTag tag) noexcept
 
EDGetToken consumes (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken consumes (TypeToGet const &id, edm::InputTag const &tag)
 
ConsumesCollector consumesCollector ()
 Use a ConsumesCollector to gather consumes information from helper functions. More...
 
template<typename ProductType , BranchType B = InEvent>
void consumesMany ()
 
void consumesMany (const TypeToGet &id)
 
template<BranchType B>
void consumesMany (const TypeToGet &id)
 
template<typename ESProduct , typename ESRecord , Transition Tr = Transition::Event>
auto esConsumes ()
 
template<typename ESProduct , typename ESRecord , Transition Tr = Transition::Event>
auto esConsumes (ESInputTag const &tag)
 
template<Transition Tr = Transition::Event>
constexpr auto esConsumes ()
 
template<Transition Tr = Transition::Event>
auto esConsumes (ESInputTag tag)
 
template<Transition Tr = Transition::Event>
ESGetTokenGeneric esConsumes (eventsetup::EventSetupRecordKey const &iRecord, eventsetup::DataKey const &iKey)
 Used with EventSetupRecord::doGet. More...
 
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > mayConsume (edm::InputTag const &tag)
 
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
void resetItemsToGetFrom (BranchType iType)
 

Detailed Description

Description: Producer of relative tracker and muon system alignment

Implementation: A sample of global muons is used for the aligning tracker and muon system relatively as "rigid bodies", i.e. determining offset and rotation (6 numbers)

Definition at line 109 of file GlobalTrackerMuonAlignment.cc.

Constructor & Destructor Documentation

◆ GlobalTrackerMuonAlignment()

GlobalTrackerMuonAlignment::GlobalTrackerMuonAlignment ( const edm::ParameterSet &  iConfig)
explicit

Definition at line 277 of file GlobalTrackerMuonAlignment.cc.

References cosmicMuonMode_, Exception, and isolatedMuonMode_.

278  : m_TkGeometryToken(esConsumes()),
279  m_MagFieldToken(esConsumes()),
280  m_globalPosToken(esConsumes()),
281  m_propToken(esConsumes(edm::ESInputTag("", iConfig.getParameter<string>("Propagator")))),
282  m_ttrhBuilderToken(esConsumes(edm::ESInputTag("", "witTrackAngle"))),
283  trackTags_(iConfig.getParameter<edm::InputTag>("tracks")),
284  muonTags_(iConfig.getParameter<edm::InputTag>("muons")),
285  gmuonTags_(iConfig.getParameter<edm::InputTag>("gmuons")),
286  smuonTags_(iConfig.getParameter<edm::InputTag>("smuons")),
287  cosmicMuonMode_(iConfig.getParameter<bool>("cosmics")),
288  isolatedMuonMode_(iConfig.getParameter<bool>("isolated")),
289  refitMuon_(iConfig.getParameter<bool>("refitmuon")),
290  refitTrack_(iConfig.getParameter<bool>("refittrack")),
291  rootOutFile_(iConfig.getUntrackedParameter<string>("rootOutFile")),
292  txtOutFile_(iConfig.getUntrackedParameter<string>("txtOutFile")),
293  writeDB_(iConfig.getUntrackedParameter<bool>("writeDB")),
294  debug_(iConfig.getUntrackedParameter<bool>("debug")),
295  defineFitter(true) {
296 #ifdef NO_FALSE_FALSE_MODE
297  if (cosmicMuonMode_ == false && isolatedMuonMode_ == false) {
298  throw cms::Exception("BadConfig") << "Muon collection not set! "
299  << "Use from GlobalTrackerMuonAlignment_test_cfg.py.";
300  }
301 #endif
302  if (cosmicMuonMode_ == true && isolatedMuonMode_ == true) {
303  throw cms::Exception("BadConfig") << "Muon collection can be either cosmic or isolated! "
304  << "Please set only one to true.";
305  }
306 }
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
GlobalTrackerMuonAlignment::m_globalPosToken
edm::ESGetToken< Alignments, GlobalPositionRcd > m_globalPosToken
Definition: GlobalTrackerMuonAlignment.cc:160
GlobalTrackerMuonAlignment::m_TkGeometryToken
edm::ESGetToken< GlobalTrackingGeometry, GlobalTrackingGeometryRecord > m_TkGeometryToken
Definition: GlobalTrackerMuonAlignment.cc:158
edm::ParameterSet::getUntrackedParameter
T getUntrackedParameter(std::string const &, T const &) const
GlobalTrackerMuonAlignment::m_ttrhBuilderToken
edm::ESGetToken< TransientTrackingRecHitBuilder, TransientRecHitRecord > m_ttrhBuilderToken
Definition: GlobalTrackerMuonAlignment.cc:162
GlobalTrackerMuonAlignment::m_MagFieldToken
edm::ESGetToken< MagneticField, IdealMagneticFieldRecord > m_MagFieldToken
Definition: GlobalTrackerMuonAlignment.cc:159
GlobalTrackerMuonAlignment::m_propToken
edm::ESGetToken< Propagator, TrackingComponentsRecord > m_propToken
Definition: GlobalTrackerMuonAlignment.cc:161

◆ ~GlobalTrackerMuonAlignment()

GlobalTrackerMuonAlignment::~GlobalTrackerMuonAlignment ( )
override

Definition at line 308 of file GlobalTrackerMuonAlignment.cc.

308  {
309  // do anything here that needs to be done at desctruction time
310  // (e.g. close files, deallocate resources etc.)
311 }

Member Function Documentation

◆ analyze()

void GlobalTrackerMuonAlignment::analyze ( const edm::Event &  iEvent,
const edm::EventSetup &  iSetup 
)
overrideprivatevirtual

Implements edm::one::EDAnalyzerBase.

Definition at line 318 of file GlobalTrackerMuonAlignment.cc.

References analyzeTrackTrajectory(), and iEvent.

318  {
319  //GlobalTrackerMuonAlignment::analyzeTrackTrack(iEvent, iSetup);
320  GlobalTrackerMuonAlignment::analyzeTrackTrajectory(iEvent, iSetup);
321 }
iEvent
int iEvent
Definition: GenABIO.cc:224
GlobalTrackerMuonAlignment::analyzeTrackTrajectory
void analyzeTrackTrajectory(const edm::Event &, const edm::EventSetup &)
Definition: GlobalTrackerMuonAlignment.cc:1254

◆ analyzeTrackTrack()

void GlobalTrackerMuonAlignment::analyzeTrackTrack ( const edm::Event &  iEvent,
const edm::EventSetup &  iSetup 
)

Definition at line 579 of file GlobalTrackerMuonAlignment.cc.

References alongMomentum, TrajectoryStateOnSurface::cartesianError(), edm::ESWatcher< T >::check(), collectionCosmic, collectionIsolated, cosmicMuonMode_, gather_cfg::cout, debug_, Vector3DBase< T, FrameTag >::dot(), MillePedeFileConverter_cfg::e, DetId::Ecal, geometryDiff::epsilon, TrajectoryStateOnSurface::freeState(), edm::EventSetup::getHandle(), TrajectoryStateOnSurface::globalMomentum(), TrajectoryStateOnSurface::globalPosition(), globalPositionRcd_, gmuonTags_, gradientGlobal(), gradientGlobalAlg(), gradientLocal(), DetId::Hcal, histo, histo101, histo102, histo11, histo12, histo13, histo14, histo15, histo16, histo17, histo18, histo19, histo2, histo20, histo21, histo22, histo23, histo24, histo25, histo26, histo27, histo28, histo29, histo3, histo30, histo31, histo32, histo33, histo34, histo35, histo4, histo5, histo6, histo7, histo8, mps_fire::i, iEvent, info(), createfilelist::int, isolatedMuonMode_, TrajectoryStateOnSurface::isValid(), iteratorEcalRcd, iteratorHcalRcd, iteratorMuonRcd, iteratorTrackerRcd, dqmiolumiharvest::j, TrajectoryStateOnSurface::localError(), TrajectoryStateOnSurface::localParameters(), visualization-live-secondInstance_cfg::m, Alignments::m_align, m_globalPosToken, m_MagFieldToken, m_propToken, m_TkGeometryToken, PV3DBase< T, PVType, FrameType >::mag(), mag(), HLT_2022v11_cff::magneticField, magneticField_, CartesianTrajectoryError::matrix(), LocalTrajectoryError::matrix(), misalignMuonL(), DetId::Muon, PDWG_BPHSkim_cff::muons, muonTags_, MuSelect, N_event, N_track, Plane::normalVector(), oppositeToMomentum, PV3DBase< T, PVType, FrameType >::perp(), GloballyPositioned< T >::position(), SmartPropagator::propagate(), TrackCandidateProducer_cfi::propagator, defaultRKPropagator::Product::propagator, DetId::rawId(), GloballyPositioned< T >::rotation(), SmartPropagator_cfi::SmartPropagator, smuonTags_, mathSSE::sqrt(), SteppingHelixPropagator_cfi::SteppingHelixPropagator, TrajectoryStateOnSurface::surface(), Surface::tangentPlane(), DetId::Tracker, trackingGeometry_, tracks, trackTags_, LocalTrajectoryParameters::vector(), watchGlobalPositionRcd_, watchMagneticFieldRecord_, watchTrackingGeometry_, x, PV3DBase< T, PVType, FrameType >::x(), TkRotation< T >::xx(), TkRotation< T >::xy(), TkRotation< T >::xz(), y, PV3DBase< T, PVType, FrameType >::y(), TkRotation< T >::yx(), TkRotation< T >::yy(), TkRotation< T >::yz(), z, PV3DBase< T, PVType, FrameType >::z(), TkRotation< T >::zx(), TkRotation< T >::zy(), and TkRotation< T >::zz().

579  {
580  using namespace edm;
581  using namespace reco;
582  //debug_ = true;
583  bool info = false;
584  bool alarm = false;
585  //bool alarm = true;
586  double PI = 3.1415927;
587 
588  cosmicMuonMode_ = true; // true: both Cosmic and IsolatedMuon are treated with 1,2 tracks
589  isolatedMuonMode_ = !cosmicMuonMode_; //true: only IsolatedMuon are treated with 1 track
590 
591  //if(N_event == 1)
592  //GlobalPositionRcdRead1::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup);
593 
594  N_event++;
595  if (info || debug_) {
596  std::cout << "----- event " << N_event << " -- tracks " << N_track << " ---";
597  if (cosmicMuonMode_)
598  std::cout << " treated as CosmicMu ";
599  if (isolatedMuonMode_)
600  std::cout << " treated as IsolatedMu ";
601  std::cout << std::endl;
602  }
603 
604  Handle<reco::TrackCollection> tracks;
605  Handle<reco::TrackCollection> muons;
606  Handle<reco::TrackCollection> gmuons;
607  Handle<reco::MuonCollection> smuons;
608 
609  if (collectionIsolated) {
610  iEvent.getByLabel("ALCARECOMuAlCalIsolatedMu", "TrackerOnly", tracks);
611  iEvent.getByLabel("ALCARECOMuAlCalIsolatedMu", "StandAlone", muons);
612  iEvent.getByLabel("ALCARECOMuAlCalIsolatedMu", "GlobalMuon", gmuons);
613  iEvent.getByLabel("ALCARECOMuAlCalIsolatedMu", "SelectedMuons", smuons);
614  } else if (collectionCosmic) {
615  iEvent.getByLabel("ALCARECOMuAlGlobalCosmics", "TrackerOnly", tracks);
616  iEvent.getByLabel("ALCARECOMuAlGlobalCosmics", "StandAlone", muons);
617  iEvent.getByLabel("ALCARECOMuAlGlobalCosmics", "GlobalMuon", gmuons);
618  iEvent.getByLabel("ALCARECOMuAlGlobalCosmics", "SelectedMuons", smuons);
619  } else {
620  iEvent.getByLabel(trackTags_, tracks);
621  iEvent.getByLabel(muonTags_, muons);
622  iEvent.getByLabel(gmuonTags_, gmuons);
623  iEvent.getByLabel(smuonTags_, smuons);
624  }
625 
626  if (debug_) {
627  std::cout << " ievBunch " << iEvent.bunchCrossing() << " runN " << (int)iEvent.run() << std::endl;
628  std::cout << " N tracks s/amu gmu selmu " << tracks->size() << " " << muons->size() << " " << gmuons->size() << " "
629  << smuons->size() << std::endl;
630  for (MuonCollection::const_iterator itMuon = smuons->begin(); itMuon != smuons->end(); ++itMuon) {
631  std::cout << " is isolatValid Matches " << itMuon->isIsolationValid() << " " << itMuon->isMatchesValid()
632  << std::endl;
633  }
634  }
635 
636  if (isolatedMuonMode_) { // ---- Only 1 Isolated Muon --------
637  if (tracks->size() != 1)
638  return;
639  if (muons->size() != 1)
640  return;
641  if (gmuons->size() != 1)
642  return;
643  if (smuons->size() != 1)
644  return;
645  }
646 
647  if (cosmicMuonMode_) { // ---- 1,2 Cosmic Muon --------
648  if (smuons->size() > 2)
649  return;
650  if (tracks->size() != smuons->size())
651  return;
652  if (muons->size() != smuons->size())
653  return;
654  if (gmuons->size() != smuons->size())
655  return;
656  }
657 
658  // ok mc_isolated_mu
659  //edm::ESHandle<TrackerGeometry> trackerGeometry;
660  //iSetup.get<TrackerDigiGeometryRecord>().get(trackerGeometry);
661  // ok mc_isolated_mu
662  //edm::ESHandle<DTGeometry> dtGeometry;
663  //iSetup.get<MuonGeometryRecord>().get(dtGeometry);
664  // don't work
665  //edm::ESHandle<CSCGeometry> cscGeometry;
666  //iSetup.get<MuonGeometryRecord>().get(cscGeometry);
667 
668  if (watchTrackingGeometry_.check(iSetup) || !trackingGeometry_) {
669  edm::ESHandle<GlobalTrackingGeometry> trackingGeometry = iSetup.getHandle(m_TkGeometryToken);
670  trackingGeometry_ = &*trackingGeometry;
671  }
672 
673  if (watchMagneticFieldRecord_.check(iSetup) || !magneticField_) {
674  edm::ESHandle<MagneticField> magneticField = iSetup.getHandle(m_MagFieldToken);
675  magneticField_ = &*magneticField;
676  }
677 
678  if (watchGlobalPositionRcd_.check(iSetup) || !globalPositionRcd_) {
679  edm::ESHandle<Alignments> globalPositionRcd = iSetup.getHandle(m_globalPosToken);
680  globalPositionRcd_ = &*globalPositionRcd;
681  for (std::vector<AlignTransform>::const_iterator i = globalPositionRcd_->m_align.begin();
682  i != globalPositionRcd_->m_align.end();
683  ++i) {
684  if (DetId(DetId::Tracker).rawId() == i->rawId())
685  iteratorTrackerRcd = i;
686  if (DetId(DetId::Muon).rawId() == i->rawId())
687  iteratorMuonRcd = i;
688  if (DetId(DetId::Ecal).rawId() == i->rawId())
689  iteratorEcalRcd = i;
690  if (DetId(DetId::Hcal).rawId() == i->rawId())
691  iteratorHcalRcd = i;
692  }
693  if (true || debug_) {
694  std::cout << "=== iteratorMuonRcd " << iteratorMuonRcd->rawId() << " ====\n"
695  << " translation " << iteratorMuonRcd->translation() << "\n"
696  << " angles " << iteratorMuonRcd->rotation().eulerAngles() << std::endl;
697  std::cout << iteratorMuonRcd->rotation() << std::endl;
698  }
699  } // end of GlobalPositionRcd
700 
701  edm::ESHandle<Propagator> propagator = iSetup.getHandle(m_propToken);
702 
703  SteppingHelixPropagator alongStHePr = SteppingHelixPropagator(magneticField_, alongMomentum);
704  SteppingHelixPropagator oppositeStHePr = SteppingHelixPropagator(magneticField_, oppositeToMomentum);
705 
706  //double tolerance = 5.e-5;
707  defaultRKPropagator::Product aprod(magneticField_, alongMomentum, 5.e-5);
708  auto& alongRKPr = aprod.propagator;
709  defaultRKPropagator::Product oprod(magneticField_, oppositeToMomentum, 5.e-5);
710  auto& oppositeRKPr = oprod.propagator;
711 
712  float epsilon = 5.;
713  SmartPropagator alongSmPr = SmartPropagator(alongRKPr, alongStHePr, magneticField_, alongMomentum, epsilon);
714  SmartPropagator oppositeSmPr =
715  SmartPropagator(oppositeRKPr, oppositeStHePr, magneticField_, oppositeToMomentum, epsilon);
716 
717  // ................................................ selected/global muon
718  //itMuon --> Jim's globalMuon
719  for (MuonCollection::const_iterator itMuon = smuons->begin(); itMuon != smuons->end(); ++itMuon) {
720  if (debug_) {
721  std::cout << " mu gM is GM Mu SaM tM " << itMuon->isGlobalMuon() << " " << itMuon->isMuon() << " "
722  << itMuon->isStandAloneMuon() << " " << itMuon->isTrackerMuon() << " " << std::endl;
723  }
724 
725  // get information about the innermost muon hit -------------------------
726  TransientTrack muTT(itMuon->outerTrack(), magneticField_, trackingGeometry_);
727  TrajectoryStateOnSurface innerMuTSOS = muTT.innermostMeasurementState();
728  TrajectoryStateOnSurface outerMuTSOS = muTT.outermostMeasurementState();
729 
730  // get information about the outermost tracker hit -----------------------
731  TransientTrack trackTT(itMuon->track(), magneticField_, trackingGeometry_);
732  TrajectoryStateOnSurface outerTrackTSOS = trackTT.outermostMeasurementState();
733  TrajectoryStateOnSurface innerTrackTSOS = trackTT.innermostMeasurementState();
734 
735  GlobalPoint pointTrackIn = innerTrackTSOS.globalPosition();
736  GlobalPoint pointTrackOut = outerTrackTSOS.globalPosition();
737  float lenghtTrack = (pointTrackOut - pointTrackIn).mag();
738  GlobalPoint pointMuonIn = innerMuTSOS.globalPosition();
739  GlobalPoint pointMuonOut = outerMuTSOS.globalPosition();
740  float lenghtMuon = (pointMuonOut - pointMuonIn).mag();
741  GlobalVector momentumTrackOut = outerTrackTSOS.globalMomentum();
742  GlobalVector momentumTrackIn = innerTrackTSOS.globalMomentum();
743  float distanceInIn = (pointTrackIn - pointMuonIn).mag();
744  float distanceInOut = (pointTrackIn - pointMuonOut).mag();
745  float distanceOutIn = (pointTrackOut - pointMuonIn).mag();
746  float distanceOutOut = (pointTrackOut - pointMuonOut).mag();
747 
748  if (debug_) {
749  std::cout << " pointTrackIn " << pointTrackIn << std::endl;
750  std::cout << " Out " << pointTrackOut << " lenght " << lenghtTrack << std::endl;
751  std::cout << " point MuonIn " << pointMuonIn << std::endl;
752  std::cout << " Out " << pointMuonOut << " lenght " << lenghtMuon << std::endl;
753  std::cout << " momeTrackIn Out " << momentumTrackIn << " " << momentumTrackOut << std::endl;
754  std::cout << " doi io ii oo " << distanceOutIn << " " << distanceInOut << " " << distanceInIn << " "
755  << distanceOutOut << std::endl;
756  }
757 
758  if (lenghtTrack < 90.)
759  continue;
760  if (lenghtMuon < 300.)
761  continue;
762  if (momentumTrackIn.mag() < 15. || momentumTrackOut.mag() < 15.)
763  continue;
764  if (trackTT.charge() != muTT.charge())
765  continue;
766 
767  if (debug_)
768  std::cout << " passed lenght momentum cuts" << std::endl;
769 
770  GlobalVector GRm, GPm, Nl, Rm, Pm, Rt, Pt, Rt0;
771  AlgebraicSymMatrix66 Cm, C0, Ce, C1;
772 
773  TrajectoryStateOnSurface extrapolationT;
774  int tsosMuonIf = 0;
775 
776  if (isolatedMuonMode_) { //------------------------------- Isolated Muon -----
777  const Surface& refSurface = innerMuTSOS.surface();
778  ConstReferenceCountingPointer<TangentPlane> tpMuLocal(refSurface.tangentPlane(innerMuTSOS.localPosition()));
779  Nl = tpMuLocal->normalVector();
780  ConstReferenceCountingPointer<TangentPlane> tpMuGlobal(refSurface.tangentPlane(innerMuTSOS.globalPosition()));
781  GlobalVector Ng = tpMuGlobal->normalVector();
782  Surface* surf = (Surface*)&refSurface;
783  const Plane* refPlane = dynamic_cast<Plane*>(surf);
784  GlobalVector Nlp = refPlane->normalVector();
785 
786  if (debug_) {
787  std::cout << " Nlocal " << Nl.x() << " " << Nl.y() << " " << Nl.z() << " "
788  << " alfa " << int(asin(Nl.x()) * 57.29578) << std::endl;
789  std::cout << " global " << Ng.x() << " " << Ng.y() << " " << Ng.z() << std::endl;
790  std::cout << " lp " << Nlp.x() << " " << Nlp.y() << " " << Nlp.z() << std::endl;
791  //std::cout<<" Nlocal Nglobal "<<Nl.x()<<" "<<Nl.y()<<" "<<Nl.z()<<" "
792  // <<Ng.x()<<" "<<Ng.y()<<" "<<Ng.z()
793  //<<" alfa "<<int(asin(Nl.x())*57.29578)<<std::endl;
794  }
795 
796  // extrapolation to innermost muon hit
797  extrapolationT = alongSmPr.propagate(outerTrackTSOS, refSurface);
798  //extrapolationT = propagator->propagate(outerTrackTSOS, refSurface);
799 
800  if (!extrapolationT.isValid()) {
801  if (false & alarm)
802  std::cout << " !!!!!!!!! Catastrophe Out-In extrapolationT.isValid "
803  //<<"\a\a\a\a\a\a\a\a"<<extrapolationT.isValid()
804  << std::endl;
805  continue;
806  }
807  tsosMuonIf = 1;
808  Rt = GlobalVector((extrapolationT.globalPosition()).x(),
809  (extrapolationT.globalPosition()).y(),
810  (extrapolationT.globalPosition()).z());
811 
812  Pt = extrapolationT.globalMomentum();
813  // global parameters of muon
814  GRm = GlobalVector(
815  (innerMuTSOS.globalPosition()).x(), (innerMuTSOS.globalPosition()).y(), (innerMuTSOS.globalPosition()).z());
816  GPm = innerMuTSOS.globalMomentum();
817  Rt0 = GlobalVector((outerTrackTSOS.globalPosition()).x(),
818  (outerTrackTSOS.globalPosition()).y(),
819  (outerTrackTSOS.globalPosition()).z());
820  Cm = AlgebraicSymMatrix66(extrapolationT.cartesianError().matrix() + innerMuTSOS.cartesianError().matrix());
821  C0 = AlgebraicSymMatrix66(outerTrackTSOS.cartesianError().matrix());
822  Ce = AlgebraicSymMatrix66(extrapolationT.cartesianError().matrix());
823  C1 = AlgebraicSymMatrix66(innerMuTSOS.cartesianError().matrix());
824 
825  } // ------------------------------- end Isolated Muon -----
826 
827  if (cosmicMuonMode_) { //------------------------------- Cosmic Muon -----
828 
829  if ((distanceOutIn <= distanceInOut) & (distanceOutIn <= distanceInIn) &
830  (distanceOutIn <= distanceOutOut)) { // ----- Out - In ------
831  if (debug_)
832  std::cout << " ----- Out - In -----" << std::endl;
833 
834  const Surface& refSurface = innerMuTSOS.surface();
835  ConstReferenceCountingPointer<TangentPlane> tpMuGlobal(refSurface.tangentPlane(innerMuTSOS.globalPosition()));
836  Nl = tpMuGlobal->normalVector();
837 
838  // extrapolation to innermost muon hit
839  extrapolationT = alongSmPr.propagate(outerTrackTSOS, refSurface);
840  //extrapolationT = propagator->propagate(outerTrackTSOS, refSurface);
841 
842  if (!extrapolationT.isValid()) {
843  if (false & alarm)
844  std::cout << " !!!!!!!!! Catastrophe Out-In extrapolationT.isValid "
845  //<<"\a\a\a\a\a\a\a\a"<<extrapolationT.isValid()
846  << std::endl;
847  continue;
848  }
849  if (debug_)
850  std::cout << " extrapolationT.isValid " << extrapolationT.isValid() << std::endl;
851 
852  tsosMuonIf = 1;
853  Rt = GlobalVector((extrapolationT.globalPosition()).x(),
854  (extrapolationT.globalPosition()).y(),
855  (extrapolationT.globalPosition()).z());
856 
857  Pt = extrapolationT.globalMomentum();
858  // global parameters of muon
859  GRm = GlobalVector(
860  (innerMuTSOS.globalPosition()).x(), (innerMuTSOS.globalPosition()).y(), (innerMuTSOS.globalPosition()).z());
861  GPm = innerMuTSOS.globalMomentum();
862  Rt0 = GlobalVector((outerTrackTSOS.globalPosition()).x(),
863  (outerTrackTSOS.globalPosition()).y(),
864  (outerTrackTSOS.globalPosition()).z());
865  Cm = AlgebraicSymMatrix66(extrapolationT.cartesianError().matrix() + innerMuTSOS.cartesianError().matrix());
866  C0 = AlgebraicSymMatrix66(outerTrackTSOS.cartesianError().matrix());
867  Ce = AlgebraicSymMatrix66(extrapolationT.cartesianError().matrix());
868  C1 = AlgebraicSymMatrix66(innerMuTSOS.cartesianError().matrix());
869 
870  if (false & debug_) {
871  //std::cout<<" ->propDir "<<propagator->propagationDirection()<<std::endl;
872  PropagationDirectionChooser Chooser;
873  std::cout << " propDirCh " << Chooser.operator()(*outerTrackTSOS.freeState(), refSurface) << " Ch == along "
874  << (alongMomentum == Chooser.operator()(*outerTrackTSOS.freeState(), refSurface)) << std::endl;
875  std::cout << " --- Nlocal " << Nl.x() << " " << Nl.y() << " " << Nl.z() << " "
876  << " alfa " << int(asin(Nl.x()) * 57.29578) << std::endl;
877  }
878  } // enf of ---- Out - In -----
879 
880  else if ((distanceInOut <= distanceInIn) & (distanceInOut <= distanceOutIn) &
881  (distanceInOut <= distanceOutOut)) { // ----- In - Out ------
882  if (debug_)
883  std::cout << " ----- In - Out -----" << std::endl;
884 
885  const Surface& refSurface = outerMuTSOS.surface();
886  ConstReferenceCountingPointer<TangentPlane> tpMuGlobal(refSurface.tangentPlane(outerMuTSOS.globalPosition()));
887  Nl = tpMuGlobal->normalVector();
888 
889  // extrapolation to outermost muon hit
890  extrapolationT = oppositeSmPr.propagate(innerTrackTSOS, refSurface);
891 
892  if (!extrapolationT.isValid()) {
893  if (false & alarm)
894  std::cout << " !!!!!!!!! Catastrophe Out-In extrapolationT.isValid "
895  << "\a\a\a\a\a\a\a\a" << extrapolationT.isValid() << std::endl;
896  continue;
897  }
898  if (debug_)
899  std::cout << " extrapolationT.isValid " << extrapolationT.isValid() << std::endl;
900 
901  tsosMuonIf = 2;
902  Rt = GlobalVector((extrapolationT.globalPosition()).x(),
903  (extrapolationT.globalPosition()).y(),
904  (extrapolationT.globalPosition()).z());
905 
906  Pt = extrapolationT.globalMomentum();
907  // global parameters of muon
908  GRm = GlobalVector(
909  (outerMuTSOS.globalPosition()).x(), (outerMuTSOS.globalPosition()).y(), (outerMuTSOS.globalPosition()).z());
910  GPm = outerMuTSOS.globalMomentum();
911  Rt0 = GlobalVector((innerTrackTSOS.globalPosition()).x(),
912  (innerTrackTSOS.globalPosition()).y(),
913  (innerTrackTSOS.globalPosition()).z());
914  Cm = AlgebraicSymMatrix66(extrapolationT.cartesianError().matrix() + outerMuTSOS.cartesianError().matrix());
915  C0 = AlgebraicSymMatrix66(innerTrackTSOS.cartesianError().matrix());
916  Ce = AlgebraicSymMatrix66(extrapolationT.cartesianError().matrix());
917  C1 = AlgebraicSymMatrix66(outerMuTSOS.cartesianError().matrix());
918 
919  if (false & debug_) {
920  //std::cout<<" ->propDir "<<propagator->propagationDirection()<<std::endl;
921  PropagationDirectionChooser Chooser;
922  std::cout << " propDirCh " << Chooser.operator()(*innerTrackTSOS.freeState(), refSurface)
923  << " Ch == oppisite "
924  << (oppositeToMomentum == Chooser.operator()(*innerTrackTSOS.freeState(), refSurface)) << std::endl;
925  std::cout << " --- Nlocal " << Nl.x() << " " << Nl.y() << " " << Nl.z() << " "
926  << " alfa " << int(asin(Nl.x()) * 57.29578) << std::endl;
927  }
928  } // enf of ---- In - Out -----
929 
930  else if ((distanceOutOut <= distanceInOut) & (distanceOutOut <= distanceInIn) &
931  (distanceOutOut <= distanceOutIn)) { // ----- Out - Out ------
932  if (debug_)
933  std::cout << " ----- Out - Out -----" << std::endl;
934 
935  // reject: momentum of track has opposite direction to muon track
936  continue;
937 
938  const Surface& refSurface = outerMuTSOS.surface();
939  ConstReferenceCountingPointer<TangentPlane> tpMuGlobal(refSurface.tangentPlane(outerMuTSOS.globalPosition()));
940  Nl = tpMuGlobal->normalVector();
941 
942  // extrapolation to outermost muon hit
943  extrapolationT = alongSmPr.propagate(outerTrackTSOS, refSurface);
944 
945  if (!extrapolationT.isValid()) {
946  if (alarm)
947  std::cout << " !!!!!!!!! Catastrophe Out-Out extrapolationT.isValid "
948  << "\a\a\a\a\a\a\a\a" << extrapolationT.isValid() << std::endl;
949  continue;
950  }
951  if (debug_)
952  std::cout << " extrapolationT.isValid " << extrapolationT.isValid() << std::endl;
953 
954  tsosMuonIf = 2;
955  Rt = GlobalVector((extrapolationT.globalPosition()).x(),
956  (extrapolationT.globalPosition()).y(),
957  (extrapolationT.globalPosition()).z());
958 
959  Pt = extrapolationT.globalMomentum();
960  // global parameters of muon
961  GRm = GlobalVector(
962  (outerMuTSOS.globalPosition()).x(), (outerMuTSOS.globalPosition()).y(), (outerMuTSOS.globalPosition()).z());
963  GPm = outerMuTSOS.globalMomentum();
964  Rt0 = GlobalVector((outerTrackTSOS.globalPosition()).x(),
965  (outerTrackTSOS.globalPosition()).y(),
966  (outerTrackTSOS.globalPosition()).z());
967  Cm = AlgebraicSymMatrix66(extrapolationT.cartesianError().matrix() + outerMuTSOS.cartesianError().matrix());
968  C0 = AlgebraicSymMatrix66(outerTrackTSOS.cartesianError().matrix());
969  Ce = AlgebraicSymMatrix66(extrapolationT.cartesianError().matrix());
970  C1 = AlgebraicSymMatrix66(outerMuTSOS.cartesianError().matrix());
971 
972  if (debug_) {
973  PropagationDirectionChooser Chooser;
974  std::cout << " propDirCh " << Chooser.operator()(*outerTrackTSOS.freeState(), refSurface) << " Ch == along "
975  << (alongMomentum == Chooser.operator()(*outerTrackTSOS.freeState(), refSurface)) << std::endl;
976  std::cout << " --- Nlocal " << Nl.x() << " " << Nl.y() << " " << Nl.z() << " "
977  << " alfa " << int(asin(Nl.x()) * 57.29578) << std::endl;
978  std::cout << " Nornal " << Nl.x() << " " << Nl.y() << " " << Nl.z() << " "
979  << " alfa " << int(asin(Nl.x()) * 57.29578) << std::endl;
980  }
981  } // enf of ---- Out - Out -----
982  else {
983  if (alarm)
984  std::cout << " ------- !!!!!!!!!!!!!!! No proper Out - In \a\a\a\a\a\a\a" << std::endl;
985  continue;
986  }
987 
988  } // ------------------------------- end Cosmic Muon -----
989 
990  if (tsosMuonIf == 0) {
991  if (info) {
992  std::cout << "No tsosMuon !!!!!!" << std::endl;
993  continue;
994  }
995  }
996  TrajectoryStateOnSurface tsosMuon;
997  if (tsosMuonIf == 1)
998  tsosMuon = muTT.innermostMeasurementState();
999  else
1000  tsosMuon = muTT.outermostMeasurementState();
1001 
1002  //GlobalTrackerMuonAlignment::misalignMuon(GRm, GPm, Nl, Rt, Rm, Pm);
1003  AlgebraicVector4 LPRm; // muon local (dx/dz, dy/dz, x, y)
1004  GlobalTrackerMuonAlignment::misalignMuonL(GRm, GPm, Nl, Rt, Rm, Pm, LPRm, extrapolationT, tsosMuon);
1005 
1006  GlobalVector resR = Rm - Rt;
1007  GlobalVector resP0 = Pm - Pt;
1008  GlobalVector resP = Pm / Pm.mag() - Pt / Pt.mag();
1009  float RelMomResidual = (Pm.mag() - Pt.mag()) / (Pt.mag() + 1.e-6);
1010  ;
1011 
1012  AlgebraicVector6 Vm;
1013  Vm(0) = resR.x();
1014  Vm(1) = resR.y();
1015  Vm(2) = resR.z();
1016  Vm(3) = resP0.x();
1017  Vm(4) = resP0.y();
1018  Vm(5) = resP0.z();
1019  float Rmuon = Rm.perp();
1020  float Zmuon = Rm.z();
1021  float alfa_x = atan2(Nl.x(), Nl.y()) * 57.29578;
1022 
1023  if (debug_) {
1024  std::cout << " Nx Ny Nz alfa_x " << Nl.x() << " " << Nl.y() << " " << Nl.z() << " " << alfa_x << std::endl;
1025  //std::cout<<" Rm "<<Rm<<std::endl<<" Rt "<<Rt<<std::endl<<" resR "<<resR<<std::endl
1026  // <<" resP "<<resP<<" dp/p "<<RelMomResidual<<std::endl;
1027  }
1028 
1029  double chi_d = 0;
1030  for (int i = 0; i <= 5; i++)
1031  chi_d += Vm(i) * Vm(i) / Cm(i, i);
1032 
1033  AlgebraicVector5 Vml(tsosMuon.localParameters().vector() - extrapolationT.localParameters().vector());
1034  AlgebraicSymMatrix55 m(tsosMuon.localError().matrix() + extrapolationT.localError().matrix());
1035  AlgebraicSymMatrix55 Cml(tsosMuon.localError().matrix() + extrapolationT.localError().matrix());
1036  bool ierrLoc = !m.Invert();
1037  if (ierrLoc && debug_ && info) {
1038  std::cout << " ==== Error inverting Local covariance matrix ==== " << std::endl;
1039  continue;
1040  }
1041  double chi_Loc = ROOT::Math::Similarity(Vml, m);
1042  if (debug_)
1043  std::cout << " chi_Loc px/pz/err " << chi_Loc << " " << Vml(1) / sqrt(Cml(1, 1)) << std::endl;
1044 
1045  if (Pt.mag() < 15.)
1046  continue;
1047  if (Pm.mag() < 5.)
1048  continue;
1049 
1050  //if(Pt.mag() < 30.) continue; // momenum cut < 30GeV
1051  //if(Pt.mag() < 60.) continue; // momenum cut < 30GeV
1052  //if(Pt.mag() > 50.) continue; // momenum cut > 50GeV
1053  //if(Pt.mag() > 100.) continue; // momenum cut > 100GeV
1054  //if(trackTT.charge() < 0) continue; // select positive charge
1055  //if(trackTT.charge() > 0) continue; // select negative charge
1056 
1057  //if(fabs(resR.x()) > 5.) continue; // strong cut X
1058  //if(fabs(resR.y()) > 5.) continue; // Y
1059  //if(fabs(resR.z()) > 5.) continue; // Z
1060  //if(fabs(resR.mag()) > 7.5) continue; // dR
1061 
1062  /*
1063  //if(fabs(RelMomResidual) > 0.5) continue;
1064  if(fabs(resR.x()) > 20.) continue;
1065  if(fabs(resR.y()) > 20.) continue;
1066  if(fabs(resR.z()) > 20.) continue;
1067  if(fabs(resR.mag()) > 30.) continue;
1068  if(fabs(resP.x()) > 0.06) continue;
1069  if(fabs(resP.y()) > 0.06) continue;
1070  if(fabs(resP.z()) > 0.06) continue;
1071  if(chi_d > 40.) continue;
1072  */
1073 
1074  // select Barrel
1075  //if(Rmuon < 400. || Rmuon > 450.) continue;
1076  //if(Zmuon < -600. || Zmuon > 600.) continue;
1077  //if(fabs(Nl.z()) > 0.95) continue;
1078  //MuSelect = " Barrel";
1079  // EndCap1
1080  //if(Rmuon < 120. || Rmuon > 450.) continue;
1081  //if(Zmuon < -720.) continue;
1082  //if(Zmuon > -580.) continue;
1083  //if(fabs(Nl.z()) < 0.95) continue;
1084  //MuSelect = " EndCap1";
1085  // EndCap2
1086  //if(Rmuon < 120. || Rmuon > 450.) continue;
1087  //if(Zmuon > 720.) continue;
1088  //if(Zmuon < 580.) continue;
1089  //if(fabs(Nl.z()) < 0.95) continue;
1090  //MuSelect = " EndCap2";
1091  // select All
1092  if (Rmuon < 120. || Rmuon > 450.)
1093  continue;
1094  if (Zmuon < -720. || Zmuon > 720.)
1095  continue;
1096  MuSelect = " Barrel+EndCaps";
1097 
1098  if (debug_)
1099  std::cout << " .............. passed all cuts" << std::endl;
1100 
1101  N_track++;
1102  // gradient and Hessian for each track
1103 
1104  GlobalTrackerMuonAlignment::gradientGlobalAlg(Rt, Pt, Rm, Nl, Cm);
1105  GlobalTrackerMuonAlignment::gradientGlobal(Rt, Pt, Rm, Pm, Nl, Cm);
1106 
1107  CLHEP::HepSymMatrix covLoc(4, 0);
1108  for (int i = 1; i <= 4; i++)
1109  for (int j = 1; j <= i; j++) {
1110  covLoc(i, j) = (tsosMuon.localError().matrix() + extrapolationT.localError().matrix())(i, j);
1111  //if(i != j) Cov(i,j) = 0.;
1112  }
1113 
1114  const Surface& refSurface = tsosMuon.surface();
1115  CLHEP::HepMatrix rotLoc(3, 3, 0);
1116  rotLoc(1, 1) = refSurface.rotation().xx();
1117  rotLoc(1, 2) = refSurface.rotation().xy();
1118  rotLoc(1, 3) = refSurface.rotation().xz();
1119 
1120  rotLoc(2, 1) = refSurface.rotation().yx();
1121  rotLoc(2, 2) = refSurface.rotation().yy();
1122  rotLoc(2, 3) = refSurface.rotation().yz();
1123 
1124  rotLoc(3, 1) = refSurface.rotation().zx();
1125  rotLoc(3, 2) = refSurface.rotation().zy();
1126  rotLoc(3, 3) = refSurface.rotation().zz();
1127 
1128  CLHEP::HepVector posLoc(3);
1129  posLoc(1) = refSurface.position().x();
1130  posLoc(2) = refSurface.position().y();
1131  posLoc(3) = refSurface.position().z();
1132 
1133  GlobalTrackerMuonAlignment::gradientLocal(Rt, Pt, Rm, Pm, Nl, covLoc, rotLoc, posLoc, LPRm);
1134 
1135  if (debug_) {
1136  std::cout << " Norm " << Nl << std::endl;
1137  std::cout << " posLoc " << posLoc.T() << std::endl;
1138  std::cout << " rotLoc " << rotLoc << std::endl;
1139  }
1140 
1141  // ----------------------------------------------------- fill histogram
1142  histo->Fill(itMuon->track()->pt());
1143 
1144  //histo2->Fill(itMuon->track()->outerP());
1145  histo2->Fill(Pt.mag());
1146  histo3->Fill((PI / 2. - itMuon->track()->outerTheta()));
1147  histo4->Fill(itMuon->track()->phi());
1148  histo5->Fill(Rmuon);
1149  histo6->Fill(Zmuon);
1150  histo7->Fill(RelMomResidual);
1151  //histo8->Fill(chi);
1152  histo8->Fill(chi_d);
1153 
1154  histo101->Fill(Zmuon, Rmuon);
1155  histo101->Fill(Rt0.z(), Rt0.perp());
1156  histo102->Fill(Rt0.x(), Rt0.y());
1157  histo102->Fill(Rm.x(), Rm.y());
1158 
1159  histo11->Fill(resR.mag());
1160  if (fabs(Nl.x()) < 0.98)
1161  histo12->Fill(resR.x());
1162  if (fabs(Nl.y()) < 0.98)
1163  histo13->Fill(resR.y());
1164  if (fabs(Nl.z()) < 0.98)
1165  histo14->Fill(resR.z());
1166  histo15->Fill(resP.x());
1167  histo16->Fill(resP.y());
1168  histo17->Fill(resP.z());
1169 
1170  if ((fabs(Nl.x()) < 0.98) && (fabs(Nl.y()) < 0.98) && (fabs(Nl.z()) < 0.98)) {
1171  histo18->Fill(sqrt(C0(0, 0)));
1172  histo19->Fill(sqrt(C1(0, 0)));
1173  histo20->Fill(sqrt(C1(0, 0) + Ce(0, 0)));
1174  }
1175  if (fabs(Nl.x()) < 0.98)
1176  histo21->Fill(Vm(0) / sqrt(Cm(0, 0)));
1177  if (fabs(Nl.y()) < 0.98)
1178  histo22->Fill(Vm(1) / sqrt(Cm(1, 1)));
1179  if (fabs(Nl.z()) < 0.98)
1180  histo23->Fill(Vm(2) / sqrt(Cm(2, 2)));
1181  histo24->Fill(Vm(3) / sqrt(C1(3, 3) + Ce(3, 3)));
1182  histo25->Fill(Vm(4) / sqrt(C1(4, 4) + Ce(4, 4)));
1183  histo26->Fill(Vm(5) / sqrt(C1(5, 5) + Ce(5, 5)));
1184  histo27->Fill(Nl.x());
1185  histo28->Fill(Nl.y());
1186  histo29->Fill(lenghtTrack);
1187  histo30->Fill(lenghtMuon);
1188  histo31->Fill(chi_Loc);
1189  histo32->Fill(Vml(1) / sqrt(Cml(1, 1)));
1190  histo33->Fill(Vml(2) / sqrt(Cml(2, 2)));
1191  histo34->Fill(Vml(3) / sqrt(Cml(3, 3)));
1192  histo35->Fill(Vml(4) / sqrt(Cml(4, 4)));
1193 
1194  if (debug_) { //--------------------------------- debug print ----------
1195 
1196  std::cout << " diag 0[ " << C0(0, 0) << " " << C0(1, 1) << " " << C0(2, 2) << " " << C0(3, 3) << " " << C0(4, 4)
1197  << " " << C0(5, 5) << " ]" << std::endl;
1198  std::cout << " diag e[ " << Ce(0, 0) << " " << Ce(1, 1) << " " << Ce(2, 2) << " " << Ce(3, 3) << " " << Ce(4, 4)
1199  << " " << Ce(5, 5) << " ]" << std::endl;
1200  std::cout << " diag 1[ " << C1(0, 0) << " " << C1(1, 1) << " " << C1(2, 2) << " " << C1(3, 3) << " " << C1(4, 4)
1201  << " " << C1(5, 5) << " ]" << std::endl;
1202  std::cout << " Rm " << Rm.x() << " " << Rm.y() << " " << Rm.z() << " Pm " << Pm.x() << " " << Pm.y() << " "
1203  << Pm.z() << std::endl;
1204  std::cout << " Rt " << Rt.x() << " " << Rt.y() << " " << Rt.z() << " Pt " << Pt.x() << " " << Pt.y() << " "
1205  << Pt.z() << std::endl;
1206  std::cout << " Nl*(Rm-Rt) " << Nl.dot(Rm - Rt) << std::endl;
1207  std::cout << " resR " << resR.x() << " " << resR.y() << " " << resR.z() << " resP " << resP.x() << " " << resP.y()
1208  << " " << resP.z() << std::endl;
1209  std::cout << " Rm-t " << (Rm - Rt).x() << " " << (Rm - Rt).y() << " " << (Rm - Rt).z() << " Pm-t "
1210  << (Pm - Pt).x() << " " << (Pm - Pt).y() << " " << (Pm - Pt).z() << std::endl;
1211  std::cout << " Vm " << Vm << std::endl;
1212  std::cout << " +- " << sqrt(Cm(0, 0)) << " " << sqrt(Cm(1, 1)) << " " << sqrt(Cm(2, 2)) << " " << sqrt(Cm(3, 3))
1213  << " " << sqrt(Cm(4, 4)) << " " << sqrt(Cm(5, 5)) << std::endl;
1214  std::cout << " Pmuon Ptrack dP/Ptrack " << itMuon->outerTrack()->p() << " " << itMuon->track()->outerP() << " "
1215  << (itMuon->outerTrack()->p() - itMuon->track()->outerP()) / itMuon->track()->outerP() << std::endl;
1216  std::cout << " cov matrix " << std::endl;
1217  std::cout << Cm << std::endl;
1218  std::cout << " diag [ " << Cm(0, 0) << " " << Cm(1, 1) << " " << Cm(2, 2) << " " << Cm(3, 3) << " " << Cm(4, 4)
1219  << " " << Cm(5, 5) << " ]" << std::endl;
1220 
1221  AlgebraicSymMatrix66 Ro;
1222  double Diag[6];
1223  for (int i = 0; i <= 5; i++)
1224  Diag[i] = sqrt(Cm(i, i));
1225  for (int i = 0; i <= 5; i++)
1226  for (int j = 0; j <= 5; j++)
1227  Ro(i, j) = Cm(i, j) / Diag[i] / Diag[j];
1228  std::cout << " correlation matrix " << std::endl;
1229  std::cout << Ro << std::endl;
1230 
1231  AlgebraicSymMatrix66 CmI;
1232  for (int i = 0; i <= 5; i++)
1233  for (int j = 0; j <= 5; j++)
1234  CmI(i, j) = Cm(i, j);
1235 
1236  bool ierr = !CmI.Invert();
1237  if (ierr) {
1238  if (alarm || debug_)
1239  std::cout << " Error inverse covariance matrix !!!!!!!!!!!" << std::endl;
1240  continue;
1241  }
1242  std::cout << " inverse cov matrix " << std::endl;
1243  std::cout << Cm << std::endl;
1244 
1245  double chi = ROOT::Math::Similarity(Vm, CmI);
1246  std::cout << " chi chi_d " << chi << " " << chi_d << std::endl;
1247  } // end of debug_ printout --------------------------------------------
1248 
1249  } // end loop on selected muons, i.e. Jim's globalMuon
1250 
1251 } //end of analyzeTrackTrack
GlobalTrackerMuonAlignment::trackingGeometry_
const GlobalTrackingGeometry * trackingGeometry_
Definition: GlobalTrackerMuonAlignment.cc:182
HLT_2022v11_cff.magneticField
magneticField
Definition: HLT_2022v11_cff.py:487
GlobalTrackerMuonAlignment::watchGlobalPositionRcd_
edm::ESWatcher< GlobalPositionRcd > watchGlobalPositionRcd_
Definition: GlobalTrackerMuonAlignment.cc:189
info
static const TGPicture * info(bool iBackgroundIsBlack)
Definition: FWCollectionSummaryWidget.cc:153
PV3DBase::perp
T perp() const
Definition: PV3DBase.h:69
mps_fire.i
i
Definition: mps_fire.py:428
SmartPropagator::propagate
TrajectoryStateOnSurface propagate(STA const &state, SUR const &surface) const
Definition: Propagator.h:50
GlobalTrackerMuonAlignment::gradientGlobalAlg
void gradientGlobalAlg(GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, AlgebraicSymMatrix66 &)
Definition: GlobalTrackerMuonAlignment.cc:2029
GlobalTrackerMuonAlignment::gradientLocal
void gradientLocal(GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, CLHEP::HepSymMatrix &, CLHEP::HepMatrix &, CLHEP::HepVector &, AlgebraicVector4 &)
Definition: GlobalTrackerMuonAlignment.cc:2305
TkRotation::xx
T xx() const
Definition: extTkRotation.h:193
TrajectoryStateOnSurface::localError
const LocalTrajectoryError & localError() const
Definition: TrajectoryStateOnSurface.h:77
Surface::tangentPlane
virtual ConstReferenceCountingPointer< TangentPlane > tangentPlane(const GlobalPoint &) const =0
PV3DBase::z
T z() const
Definition: PV3DBase.h:61
TrackCandidateProducer_cfi.propagator
propagator
Definition: TrackCandidateProducer_cfi.py:17
Vector3DBase::dot
PreciseFloatType< T, U >::Type dot(const Vector3DBase< U, FrameTag > &v) const
Definition: Vector3DBase.h:99
TkRotation::xy
T xy() const
Definition: extTkRotation.h:194
GlobalTrackerMuonAlignment::m_globalPosToken
edm::ESGetToken< Alignments, GlobalPositionRcd > m_globalPosToken
Definition: GlobalTrackerMuonAlignment.cc:160
GlobalTrackerMuonAlignment::m_TkGeometryToken
edm::ESGetToken< GlobalTrackingGeometry, GlobalTrackingGeometryRecord > m_TkGeometryToken
Definition: GlobalTrackerMuonAlignment.cc:158
SteppingHelixPropagator_cfi.SteppingHelixPropagator
SteppingHelixPropagator
Definition: SteppingHelixPropagator_cfi.py:3
TkRotation::zz
T zz() const
Definition: extTkRotation.h:201
GlobalTrackerMuonAlignment::iteratorHcalRcd
std::vector< AlignTransform >::const_iterator iteratorHcalRcd
Definition: GlobalTrackerMuonAlignment.cc:217
TkRotation::yy
T yy() const
Definition: extTkRotation.h:197
TrajectoryStateOnSurface::localParameters
const LocalTrajectoryParameters & localParameters() const
Definition: TrajectoryStateOnSurface.h:73
TkRotation::yz
T yz() const
Definition: extTkRotation.h:198
MillePedeFileConverter_cfg.e
e
Definition: MillePedeFileConverter_cfg.py:37
TrajectoryStateOnSurface::surface
const SurfaceType & surface() const
Definition: TrajectoryStateOnSurface.h:78
GlobalTrackerMuonAlignment::gradientGlobal
void gradientGlobal(GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, AlgebraicSymMatrix66 &)
Definition: GlobalTrackerMuonAlignment.cc:2120
Alignments::m_align
std::vector< AlignTransform > m_align
Definition: Alignments.h:19
Plane
Definition: Plane.h:16
TrajectoryStateOnSurface::cartesianError
const CartesianTrajectoryError cartesianError() const
Definition: TrajectoryStateOnSurface.h:71
TkRotation::zx
T zx() const
Definition: extTkRotation.h:199
PV3DBase::x
T x() const
Definition: PV3DBase.h:59
PV3DBase::y
T y() const
Definition: PV3DBase.h:60
iEvent
int iEvent
Definition: GenABIO.cc:224
TrajectoryStateOnSurface::globalPosition
GlobalPoint globalPosition() const
Definition: TrajectoryStateOnSurface.h:65
GlobalTrackerMuonAlignment::watchTrackingGeometry_
edm::ESWatcher< GlobalTrackingGeometryRecord > watchTrackingGeometry_
Definition: GlobalTrackerMuonAlignment.cc:180
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
LocalTrajectoryParameters::vector
AlgebraicVector5 vector() const
Definition: LocalTrajectoryParameters.h:120
visualization-live-secondInstance_cfg.m
m
Definition: visualization-live-secondInstance_cfg.py:80
TkRotation::zy
T zy() const
Definition: extTkRotation.h:200
PV3DBase::mag
T mag() const
Definition: PV3DBase.h:64
GlobalTrackerMuonAlignment::iteratorMuonRcd
std::vector< AlignTransform >::const_iterator iteratorMuonRcd
Definition: GlobalTrackerMuonAlignment.cc:215
AlgebraicVector5
ROOT::Math::SVector< double, 5 > AlgebraicVector5
Definition: AlgebraicROOTObjects.h:14
AlgebraicVector4
ROOT::Math::SVector< double, 4 > AlgebraicVector4
Definition: AlgebraicROOTObjects.h:13
edm::EventSetup::getHandle
ESHandle< T > getHandle(const ESGetToken< T, R > &iToken) const
Definition: EventSetup.h:151
GlobalTrackerMuonAlignment::misalignMuonL
void misalignMuonL(GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, AlgebraicVector4 &, TrajectoryStateOnSurface &, TrajectoryStateOnSurface &)
Definition: GlobalTrackerMuonAlignment.cc:2712
AlgebraicSymMatrix66
ROOT::Math::SMatrix< double, 6, 6, ROOT::Math::MatRepSym< double, 6 > > AlgebraicSymMatrix66
Definition: AlgebraicROOTObjects.h:24
GlobalTrackerMuonAlignment::iteratorEcalRcd
std::vector< AlignTransform >::const_iterator iteratorEcalRcd
Definition: GlobalTrackerMuonAlignment.cc:216
DetId
Definition: DetId.h:17
mag
T mag() const
The vector magnitude. Equivalent to sqrt(vec.mag2())
Definition: Basic3DVectorLD.h:127
tracks
auto const & tracks
cannot be loose
Definition: CAHitNtupletGeneratorKernelsImpl.h:148
GloballyPositioned::position
const PositionType & position() const
Definition: GloballyPositioned.h:36
DetId::rawId
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:57
SmartPropagator_cfi.SmartPropagator
SmartPropagator
Definition: SmartPropagator_cfi.py:3
AlgebraicSymMatrix55
ROOT::Math::SMatrix< double, 5, 5, ROOT::Math::MatRepSym< double, 5 > > AlgebraicSymMatrix55
Definition: AlgebraicROOTObjects.h:23
GlobalTrackerMuonAlignment::watchMagneticFieldRecord_
edm::ESWatcher< IdealMagneticFieldRecord > watchMagneticFieldRecord_
Definition: GlobalTrackerMuonAlignment.cc:184
TrajectoryStateOnSurface::globalMomentum
GlobalVector globalMomentum() const
Definition: TrajectoryStateOnSurface.h:66
edm::ESWatcher::check
bool check(const edm::EventSetup &iSetup)
Definition: ESWatcher.h:57
GlobalTrackerMuonAlignment::m_MagFieldToken
edm::ESGetToken< MagneticField, IdealMagneticFieldRecord > m_MagFieldToken
Definition: GlobalTrackerMuonAlignment.cc:159
reco
fixed size matrix
Definition: AlignmentAlgorithmBase.h:46
edm
HLT enums.
Definition: AlignableModifier.h:19
LocalTrajectoryError::matrix
const AlgebraicSymMatrix55 & matrix() const
Definition: LocalTrajectoryError.h:60
CartesianTrajectoryError::matrix
const AlgebraicSymMatrix66 & matrix() const
Definition: CartesianTrajectoryError.h:28
TrajectoryStateOnSurface::freeState
FreeTrajectoryState const * freeState(bool withErrors=true) const
Definition: TrajectoryStateOnSurface.h:58
AlgebraicVector6
ROOT::Math::SVector< double, 6 > AlgebraicVector6
Definition: AlgebraicROOTObjects.h:15
Plane::normalVector
GlobalVector normalVector() const
Definition: Plane.h:41
GlobalTrackerMuonAlignment::iteratorTrackerRcd
std::vector< AlignTransform >::const_iterator iteratorTrackerRcd
Definition: GlobalTrackerMuonAlignment.cc:214
GloballyPositioned::rotation
const RotationType & rotation() const
Definition: GloballyPositioned.h:38
PDWG_BPHSkim_cff.muons
muons
Definition: PDWG_BPHSkim_cff.py:47
gather_cfg.cout
cout
Definition: gather_cfg.py:144
GlobalTrackerMuonAlignment::m_propToken
edm::ESGetToken< Propagator, TrackingComponentsRecord > m_propToken
Definition: GlobalTrackerMuonAlignment.cc:161
TkRotation::xz
T xz() const
Definition: extTkRotation.h:195
TkRotation::yx
T yx() const
Definition: extTkRotation.h:196
GlobalVector
Global3DVector GlobalVector
Definition: GlobalVector.h:10

◆ analyzeTrackTrajectory()

void GlobalTrackerMuonAlignment::analyzeTrackTrajectory ( const edm::Event &  iEvent,
const edm::EventSetup &  iSetup 
)

Definition at line 1254 of file GlobalTrackerMuonAlignment.cc.

References alongMomentum, TrajectoryStateOnSurface::cartesianError(), edm::ESWatcher< T >::check(), Chi2MeasurementEstimator_cfi::Chi2MeasurementEstimator, collectionCosmic, collectionIsolated, cosmicMuonMode_, gather_cfg::cout, debug_, debugTrajectorySOS(), defineFitter, Vector3DBase< T, FrameTag >::dot(), MillePedeFileConverter_cfg::e, DetId::Ecal, geometryDiff::epsilon, TrajectoryStateOnSurface::freeState(), edm::EventSetup::getHandle(), TrajectoryStateOnSurface::globalMomentum(), TrajectoryStateOnSurface::globalPosition(), globalPositionRcd_, gmuonTags_, gradientGlobal(), gradientGlobalAlg(), gradientLocal(), DetId::Hcal, histo, histo101, histo102, histo11, histo12, histo13, histo14, histo15, histo16, histo17, histo18, histo19, histo2, histo20, histo21, histo22, histo23, histo24, histo25, histo26, histo27, histo28, histo29, histo3, histo30, histo31, histo32, histo33, histo34, histo35, histo4, histo5, histo6, histo7, histo8, mps_fire::i, iEvent, info(), createfilelist::int, isolatedMuonMode_, edm::ESHandleBase::isValid(), TrajectoryStateOnSurface::isValid(), iteratorEcalRcd, iteratorHcalRcd, iteratorMuonRcd, iteratorTrackerRcd, dqmiolumiharvest::j, TrajectoryStateOnSurface::localError(), TrajectoryStateOnSurface::localParameters(), visualization-live-secondInstance_cfg::m, Alignments::m_align, m_globalPosToken, m_MagFieldToken, m_propToken, m_TkGeometryToken, m_ttrhBuilderToken, PV3DBase< T, PVType, FrameType >::mag(), mag(), HLT_2022v11_cff::magneticField, magneticField_, CartesianTrajectoryError::matrix(), LocalTrajectoryError::matrix(), misalignMuonL(), DetId::Muon, muonFitter(), PDWG_BPHSkim_cff::muons, muonTags_, MuRHBuilder, MuSelect, N_event, N_track, Plane::normalVector(), oppositeToMomentum, PV3DBase< T, PVType, FrameType >::perp(), GloballyPositioned< T >::position(), SmartPropagator::propagate(), TrackCandidateProducer_cfi::propagator, defaultRKPropagator::Product::propagator, DetId::rawId(), refitMuon_, refitTrack_, GloballyPositioned< T >::rotation(), SmartPropagator_cfi::SmartPropagator, smuonTags_, mathSSE::sqrt(), SteppingHelixPropagator_cfi::SteppingHelixPropagator, TrajectoryStateOnSurface::surface(), Surface::tangentPlane(), theFitter, theFitterOp, theSmoother, theSmootherOp, theTrackingGeometry, DetId::Tracker, trackFitter(), trackingGeometry_, tracks, trackTags_, TTRHBuilder, LocalTrajectoryParameters::vector(), watchGlobalPositionRcd_, watchMagneticFieldRecord_, watchTrackingGeometry_, x, PV3DBase< T, PVType, FrameType >::x(), TkRotation< T >::xx(), TkRotation< T >::xy(), TkRotation< T >::xz(), y, PV3DBase< T, PVType, FrameType >::y(), TkRotation< T >::yx(), TkRotation< T >::yy(), TkRotation< T >::yz(), z, PV3DBase< T, PVType, FrameType >::z(), TkRotation< T >::zx(), TkRotation< T >::zy(), and TkRotation< T >::zz().

Referenced by analyze().

1254  {
1255  using namespace edm;
1256  using namespace reco;
1257  //debug_ = true;
1258  bool info = false;
1259  bool alarm = false;
1260  //bool alarm = true;
1261  double PI = 3.1415927;
1262 
1263  //-*-*-*-*-*-*-
1264  cosmicMuonMode_ = true; // true: both Cosmic and IsolatedMuon are treated with 1,2 tracks
1265  //cosmicMuonMode_ = false; // true: both Cosmic and IsolatedMuon are treated with 1,2 tracks
1266  //-*-*-*-*-*-*-
1267  isolatedMuonMode_ = !cosmicMuonMode_; //true: only IsolatedMuon are treated with 1 track
1268 
1269  N_event++;
1270  if (info || debug_) {
1271  std::cout << "----- event " << N_event << " -- tracks " << N_track << " ---";
1272  if (cosmicMuonMode_)
1273  std::cout << " treated as CosmicMu ";
1274  if (isolatedMuonMode_)
1275  std::cout << " treated as IsolatedMu ";
1276  std::cout << std::endl;
1277  }
1278 
1279  Handle<reco::TrackCollection> tracks;
1280  Handle<reco::TrackCollection> muons;
1281  Handle<reco::TrackCollection> gmuons;
1282  Handle<reco::MuonCollection> smuons;
1283 
1284  if (collectionIsolated) {
1285  iEvent.getByLabel("ALCARECOMuAlCalIsolatedMu", "TrackerOnly", tracks);
1286  iEvent.getByLabel("ALCARECOMuAlCalIsolatedMu", "StandAlone", muons);
1287  iEvent.getByLabel("ALCARECOMuAlCalIsolatedMu", "GlobalMuon", gmuons);
1288  iEvent.getByLabel("ALCARECOMuAlCalIsolatedMu", "SelectedMuons", smuons);
1289  } else if (collectionCosmic) {
1290  iEvent.getByLabel("ALCARECOMuAlGlobalCosmics", "TrackerOnly", tracks);
1291  iEvent.getByLabel("ALCARECOMuAlGlobalCosmics", "StandAlone", muons);
1292  iEvent.getByLabel("ALCARECOMuAlGlobalCosmics", "GlobalMuon", gmuons);
1293  iEvent.getByLabel("ALCARECOMuAlGlobalCosmics", "SelectedMuons", smuons);
1294  } else {
1295  iEvent.getByLabel(trackTags_, tracks);
1296  iEvent.getByLabel(muonTags_, muons);
1297  iEvent.getByLabel(gmuonTags_, gmuons);
1298  iEvent.getByLabel(smuonTags_, smuons);
1299  }
1300 
1301  if (debug_) {
1302  std::cout << " ievBunch " << iEvent.bunchCrossing() << " runN " << (int)iEvent.run() << std::endl;
1303  std::cout << " N tracks s/amu gmu selmu " << tracks->size() << " " << muons->size() << " " << gmuons->size() << " "
1304  << smuons->size() << std::endl;
1305  for (MuonCollection::const_iterator itMuon = smuons->begin(); itMuon != smuons->end(); ++itMuon) {
1306  std::cout << " is isolatValid Matches " << itMuon->isIsolationValid() << " " << itMuon->isMatchesValid()
1307  << std::endl;
1308  }
1309  }
1310 
1311  if (isolatedMuonMode_) { // ---- Only 1 Isolated Muon --------
1312  if (tracks->size() != 1)
1313  return;
1314  if (muons->size() != 1)
1315  return;
1316  if (gmuons->size() != 1)
1317  return;
1318  if (smuons->size() != 1)
1319  return;
1320  }
1321 
1322  if (cosmicMuonMode_) { // ---- 1,2 Cosmic Muon --------
1323  if (smuons->size() > 2)
1324  return;
1325  if (tracks->size() != smuons->size())
1326  return;
1327  if (muons->size() != smuons->size())
1328  return;
1329  if (gmuons->size() != smuons->size())
1330  return;
1331  }
1332 
1333  // ok mc_isolated_mu
1334  //edm::ESHandle<TrackerGeometry> trackerGeometry;
1335  //iSetup.get<TrackerDigiGeometryRecord>().get(trackerGeometry);
1336  // ok mc_isolated_mu
1337  //edm::ESHandle<DTGeometry> dtGeometry;
1338  //iSetup.get<MuonGeometryRecord>().get(dtGeometry);
1339  // don't work
1340  //edm::ESHandle<CSCGeometry> cscGeometry;
1341  //iSetup.get<MuonGeometryRecord>().get(cscGeometry);
1342 
1343  if (watchTrackingGeometry_.check(iSetup) || !trackingGeometry_) {
1344  edm::ESHandle<GlobalTrackingGeometry> trackingGeometry = iSetup.getHandle(m_TkGeometryToken);
1345  trackingGeometry_ = &*trackingGeometry;
1346  theTrackingGeometry = trackingGeometry;
1347  }
1348 
1349  if (watchMagneticFieldRecord_.check(iSetup) || !magneticField_) {
1350  edm::ESHandle<MagneticField> magneticField = iSetup.getHandle(m_MagFieldToken);
1351  magneticField_ = &*magneticField;
1352  }
1353 
1354  if (watchGlobalPositionRcd_.check(iSetup) || !globalPositionRcd_) {
1355  edm::ESHandle<Alignments> globalPositionRcd = iSetup.getHandle(m_globalPosToken);
1356  globalPositionRcd_ = &*globalPositionRcd;
1357  for (std::vector<AlignTransform>::const_iterator i = globalPositionRcd_->m_align.begin();
1358  i != globalPositionRcd_->m_align.end();
1359  ++i) {
1360  if (DetId(DetId::Tracker).rawId() == i->rawId())
1361  iteratorTrackerRcd = i;
1362  if (DetId(DetId::Muon).rawId() == i->rawId())
1363  iteratorMuonRcd = i;
1364  if (DetId(DetId::Ecal).rawId() == i->rawId())
1365  iteratorEcalRcd = i;
1366  if (DetId(DetId::Hcal).rawId() == i->rawId())
1367  iteratorHcalRcd = i;
1368  }
1369  if (debug_) {
1370  std::cout << "=== iteratorTrackerRcd " << iteratorTrackerRcd->rawId() << " ====\n"
1371  << " translation " << iteratorTrackerRcd->translation() << "\n"
1372  << " angles " << iteratorTrackerRcd->rotation().eulerAngles() << std::endl;
1373  std::cout << iteratorTrackerRcd->rotation() << std::endl;
1374  std::cout << "=== iteratorMuonRcd " << iteratorMuonRcd->rawId() << " ====\n"
1375  << " translation " << iteratorMuonRcd->translation() << "\n"
1376  << " angles " << iteratorMuonRcd->rotation().eulerAngles() << std::endl;
1377  std::cout << iteratorMuonRcd->rotation() << std::endl;
1378  }
1379  } // end of GlobalPositionRcd
1380 
1381  edm::ESHandle<Propagator> propagator = iSetup.getHandle(m_propToken);
1382 
1383  SteppingHelixPropagator alongStHePr = SteppingHelixPropagator(magneticField_, alongMomentum);
1384  SteppingHelixPropagator oppositeStHePr = SteppingHelixPropagator(magneticField_, oppositeToMomentum);
1385 
1386  //double tolerance = 5.e-5;
1387  defaultRKPropagator::Product aprod(magneticField_, alongMomentum, 5.e-5);
1388  auto& alongRKPr = aprod.propagator;
1389  defaultRKPropagator::Product oprod(magneticField_, oppositeToMomentum, 5.e-5);
1390  auto& oppositeRKPr = oprod.propagator;
1391 
1392  float epsilon = 5.;
1393  SmartPropagator alongSmPr = SmartPropagator(alongRKPr, alongStHePr, magneticField_, alongMomentum, epsilon);
1394  SmartPropagator oppositeSmPr =
1395  SmartPropagator(oppositeRKPr, oppositeStHePr, magneticField_, oppositeToMomentum, epsilon);
1396 
1397  if (defineFitter) {
1398  if (debug_)
1399  std::cout << " ............... DEFINE FITTER ..................." << std::endl;
1400  KFUpdator* theUpdator = new KFUpdator();
1401  //Chi2MeasurementEstimator* theEstimator = new Chi2MeasurementEstimator(30);
1402  Chi2MeasurementEstimator* theEstimator = new Chi2MeasurementEstimator(100000, 100000);
1403  theFitter = new KFTrajectoryFitter(alongSmPr, *theUpdator, *theEstimator);
1404  theSmoother = new KFTrajectorySmoother(alongSmPr, *theUpdator, *theEstimator);
1405  theFitterOp = new KFTrajectoryFitter(oppositeSmPr, *theUpdator, *theEstimator);
1406  theSmootherOp = new KFTrajectorySmoother(oppositeSmPr, *theUpdator, *theEstimator);
1407 
1408  edm::ESHandle<TransientTrackingRecHitBuilder> builder = iSetup.getHandle(m_ttrhBuilderToken);
1409  this->TTRHBuilder = &(*builder);
1410  this->MuRHBuilder = new MuonTransientTrackingRecHitBuilder(theTrackingGeometry);
1411  if (debug_) {
1412  std::cout << " theTrackingGeometry.isValid() " << theTrackingGeometry.isValid() << std::endl;
1413  std::cout << "get also the MuonTransientTrackingRecHitBuilder"
1414  << "\n";
1415  std::cout << "get also the TransientTrackingRecHitBuilder"
1416  << "\n";
1417  }
1418  defineFitter = false;
1419  }
1420 
1421  // ................................................ selected/global muon
1422  //itMuon --> Jim's globalMuon
1423  for (MuonCollection::const_iterator itMuon = smuons->begin(); itMuon != smuons->end(); ++itMuon) {
1424  if (debug_) {
1425  std::cout << " mu gM is GM Mu SaM tM " << itMuon->isGlobalMuon() << " " << itMuon->isMuon() << " "
1426  << itMuon->isStandAloneMuon() << " " << itMuon->isTrackerMuon() << " " << std::endl;
1427  }
1428 
1429  // get information about the innermost muon hit -------------------------
1430  TransientTrack muTT(itMuon->outerTrack(), magneticField_, trackingGeometry_);
1431  TrajectoryStateOnSurface innerMuTSOS = muTT.innermostMeasurementState();
1432  TrajectoryStateOnSurface outerMuTSOS = muTT.outermostMeasurementState();
1433 
1434  // get information about the outermost tracker hit -----------------------
1435  TransientTrack trackTT(itMuon->track(), magneticField_, trackingGeometry_);
1436  TrajectoryStateOnSurface outerTrackTSOS = trackTT.outermostMeasurementState();
1437  TrajectoryStateOnSurface innerTrackTSOS = trackTT.innermostMeasurementState();
1438 
1439  GlobalPoint pointTrackIn = innerTrackTSOS.globalPosition();
1440  GlobalPoint pointTrackOut = outerTrackTSOS.globalPosition();
1441  float lenghtTrack = (pointTrackOut - pointTrackIn).mag();
1442  GlobalPoint pointMuonIn = innerMuTSOS.globalPosition();
1443  GlobalPoint pointMuonOut = outerMuTSOS.globalPosition();
1444  float lenghtMuon = (pointMuonOut - pointMuonIn).mag();
1445  GlobalVector momentumTrackOut = outerTrackTSOS.globalMomentum();
1446  GlobalVector momentumTrackIn = innerTrackTSOS.globalMomentum();
1447  float distanceInIn = (pointTrackIn - pointMuonIn).mag();
1448  float distanceInOut = (pointTrackIn - pointMuonOut).mag();
1449  float distanceOutIn = (pointTrackOut - pointMuonIn).mag();
1450  float distanceOutOut = (pointTrackOut - pointMuonOut).mag();
1451 
1452  if (debug_) {
1453  std::cout << " pointTrackIn " << pointTrackIn << std::endl;
1454  std::cout << " Out " << pointTrackOut << " lenght " << lenghtTrack << std::endl;
1455  std::cout << " point MuonIn " << pointMuonIn << std::endl;
1456  std::cout << " Out " << pointMuonOut << " lenght " << lenghtMuon << std::endl;
1457  std::cout << " momeTrackIn Out " << momentumTrackIn << " " << momentumTrackOut << std::endl;
1458  std::cout << " doi io ii oo " << distanceOutIn << " " << distanceInOut << " " << distanceInIn << " "
1459  << distanceOutOut << std::endl;
1460  }
1461 
1462  if (lenghtTrack < 90.)
1463  continue;
1464  if (lenghtMuon < 300.)
1465  continue;
1466  if (innerMuTSOS.globalMomentum().mag() < 5. || outerMuTSOS.globalMomentum().mag() < 5.)
1467  continue;
1468  if (momentumTrackIn.mag() < 15. || momentumTrackOut.mag() < 15.)
1469  continue;
1470  if (trackTT.charge() != muTT.charge())
1471  continue;
1472 
1473  if (debug_)
1474  std::cout << " passed lenght momentum cuts" << std::endl;
1475 
1476  GlobalVector GRm, GPm, Nl, Rm, Pm, Rt, Pt, Rt0;
1477  AlgebraicSymMatrix66 Cm, C0, Ce, C1;
1478 
1479  //extrapolationT = propagator->propagate(outerTrackTSOS, refSurface);
1480  TrajectoryStateOnSurface extrapolationT;
1481  int tsosMuonIf = 0;
1482 
1483  TrajectoryStateOnSurface muonFittedTSOS;
1484  TrajectoryStateOnSurface trackFittedTSOS;
1485 
1486  if (isolatedMuonMode_) { //------------------------------- Isolated Muon --- Out-In --
1487  if (debug_)
1488  std::cout << " ------ Isolated (out-in) ----- " << std::endl;
1489  const Surface& refSurface = innerMuTSOS.surface();
1490  ConstReferenceCountingPointer<TangentPlane> tpMuLocal(refSurface.tangentPlane(innerMuTSOS.localPosition()));
1491  Nl = tpMuLocal->normalVector();
1492  ConstReferenceCountingPointer<TangentPlane> tpMuGlobal(refSurface.tangentPlane(innerMuTSOS.globalPosition()));
1493  GlobalVector Ng = tpMuGlobal->normalVector();
1494  Surface* surf = (Surface*)&refSurface;
1495  const Plane* refPlane = dynamic_cast<Plane*>(surf);
1496  GlobalVector Nlp = refPlane->normalVector();
1497 
1498  if (debug_) {
1499  std::cout << " Nlocal " << Nl.x() << " " << Nl.y() << " " << Nl.z() << " "
1500  << " alfa " << int(asin(Nl.x()) * 57.29578) << std::endl;
1501  std::cout << " global " << Ng.x() << " " << Ng.y() << " " << Ng.z() << std::endl;
1502  std::cout << " lp " << Nlp.x() << " " << Nlp.y() << " " << Nlp.z() << std::endl;
1503  //std::cout<<" Nlocal Nglobal "<<Nl.x()<<" "<<Nl.y()<<" "<<Nl.z()<<" "
1504  // <<Ng.x()<<" "<<Ng.y()<<" "<<Ng.z()
1505  //<<" alfa "<<int(asin(Nl.x())*57.29578)<<std::endl;
1506  }
1507 
1508  // extrapolation to innermost muon hit
1509 
1510  //extrapolationT = alongSmPr.propagate(outerTrackTSOS, refSurface);
1511  if (!refitTrack_)
1512  extrapolationT = alongSmPr.propagate(outerTrackTSOS, refSurface);
1513  else {
1514  GlobalTrackerMuonAlignment::trackFitter(itMuon->track(), trackTT, alongMomentum, trackFittedTSOS);
1515  if (trackFittedTSOS.isValid())
1516  extrapolationT = alongSmPr.propagate(trackFittedTSOS, refSurface);
1517  }
1518 
1519  if (!extrapolationT.isValid()) {
1520  if (false & alarm)
1521  std::cout << " !!!!!!!!! Catastrophe Out-In extrapolationT.isValid "
1522  //<<"\a\a\a\a\a\a\a\a"<<extrapolationT.isValid()
1523  << std::endl;
1524  continue;
1525  }
1526  tsosMuonIf = 1;
1527  Rt = GlobalVector((extrapolationT.globalPosition()).x(),
1528  (extrapolationT.globalPosition()).y(),
1529  (extrapolationT.globalPosition()).z());
1530 
1531  Pt = extrapolationT.globalMomentum();
1532  // global parameters of muon
1533  GRm = GlobalVector(
1534  (innerMuTSOS.globalPosition()).x(), (innerMuTSOS.globalPosition()).y(), (innerMuTSOS.globalPosition()).z());
1535  GPm = innerMuTSOS.globalMomentum();
1536 
1537  Rt0 = GlobalVector((outerTrackTSOS.globalPosition()).x(),
1538  (outerTrackTSOS.globalPosition()).y(),
1539  (outerTrackTSOS.globalPosition()).z());
1540  Cm = AlgebraicSymMatrix66(extrapolationT.cartesianError().matrix() + innerMuTSOS.cartesianError().matrix());
1541  C0 = AlgebraicSymMatrix66(outerTrackTSOS.cartesianError().matrix());
1542  Ce = AlgebraicSymMatrix66(extrapolationT.cartesianError().matrix());
1543  C1 = AlgebraicSymMatrix66(innerMuTSOS.cartesianError().matrix());
1544 
1545  if (refitMuon_)
1546  GlobalTrackerMuonAlignment::muonFitter(itMuon->outerTrack(), muTT, oppositeToMomentum, muonFittedTSOS);
1547 
1548  } // ------------------------------- end Isolated Muon -- Out - In ---
1549 
1550  if (cosmicMuonMode_) { //------------------------------- Cosmic Muon -----
1551 
1552  if ((distanceOutIn <= distanceInOut) & (distanceOutIn <= distanceInIn) &
1553  (distanceOutIn <= distanceOutOut)) { // ----- Out - In ------
1554  if (debug_)
1555  std::cout << " ----- Out - In -----" << std::endl;
1556 
1557  const Surface& refSurface = innerMuTSOS.surface();
1558  ConstReferenceCountingPointer<TangentPlane> tpMuGlobal(refSurface.tangentPlane(innerMuTSOS.globalPosition()));
1559  Nl = tpMuGlobal->normalVector();
1560 
1561  // extrapolation to innermost muon hit
1562  //extrapolationT = alongSmPr.propagate(outerTrackTSOS, refSurface);
1563  if (!refitTrack_)
1564  extrapolationT = alongSmPr.propagate(outerTrackTSOS, refSurface);
1565  else {
1566  GlobalTrackerMuonAlignment::trackFitter(itMuon->track(), trackTT, alongMomentum, trackFittedTSOS);
1567  if (trackFittedTSOS.isValid())
1568  extrapolationT = alongSmPr.propagate(trackFittedTSOS, refSurface);
1569  }
1570 
1571  if (!extrapolationT.isValid()) {
1572  if (false & alarm)
1573  std::cout << " !!!!!!!!! Catastrophe Out-In extrapolationT.isValid "
1574  //<<"\a\a\a\a\a\a\a\a"<<extrapolationT.isValid()
1575  << std::endl;
1576  continue;
1577  }
1578  if (debug_)
1579  std::cout << " extrapolationT.isValid " << extrapolationT.isValid() << std::endl;
1580 
1581  tsosMuonIf = 1;
1582  Rt = GlobalVector((extrapolationT.globalPosition()).x(),
1583  (extrapolationT.globalPosition()).y(),
1584  (extrapolationT.globalPosition()).z());
1585 
1586  Pt = extrapolationT.globalMomentum();
1587  // global parameters of muon
1588  GRm = GlobalVector(
1589  (innerMuTSOS.globalPosition()).x(), (innerMuTSOS.globalPosition()).y(), (innerMuTSOS.globalPosition()).z());
1590  GPm = innerMuTSOS.globalMomentum();
1591  Rt0 = GlobalVector((outerTrackTSOS.globalPosition()).x(),
1592  (outerTrackTSOS.globalPosition()).y(),
1593  (outerTrackTSOS.globalPosition()).z());
1594  Cm = AlgebraicSymMatrix66(extrapolationT.cartesianError().matrix() + innerMuTSOS.cartesianError().matrix());
1595  C0 = AlgebraicSymMatrix66(outerTrackTSOS.cartesianError().matrix());
1596  Ce = AlgebraicSymMatrix66(extrapolationT.cartesianError().matrix());
1597  C1 = AlgebraicSymMatrix66(innerMuTSOS.cartesianError().matrix());
1598 
1599  if (refitMuon_)
1600  GlobalTrackerMuonAlignment::muonFitter(itMuon->outerTrack(), muTT, oppositeToMomentum, muonFittedTSOS);
1601 
1602  if (false & debug_) {
1603  //std::cout<<" ->propDir "<<propagator->propagationDirection()<<std::endl;
1604  PropagationDirectionChooser Chooser;
1605  std::cout << " propDirCh " << Chooser.operator()(*outerTrackTSOS.freeState(), refSurface) << " Ch == along "
1606  << (alongMomentum == Chooser.operator()(*outerTrackTSOS.freeState(), refSurface)) << std::endl;
1607  std::cout << " --- Nlocal " << Nl.x() << " " << Nl.y() << " " << Nl.z() << " "
1608  << " alfa " << int(asin(Nl.x()) * 57.29578) << std::endl;
1609  }
1610  } // enf of ---- Out - In -----
1611 
1612  else if ((distanceInOut <= distanceInIn) & (distanceInOut <= distanceOutIn) &
1613  (distanceInOut <= distanceOutOut)) { // ----- In - Out ------
1614  if (debug_)
1615  std::cout << " ----- In - Out -----" << std::endl;
1616 
1617  const Surface& refSurface = outerMuTSOS.surface();
1618  ConstReferenceCountingPointer<TangentPlane> tpMuGlobal(refSurface.tangentPlane(outerMuTSOS.globalPosition()));
1619  Nl = tpMuGlobal->normalVector();
1620 
1621  // extrapolation to outermost muon hit
1622  //extrapolationT = oppositeSmPr.propagate(innerTrackTSOS, refSurface);
1623  if (!refitTrack_)
1624  extrapolationT = oppositeSmPr.propagate(innerTrackTSOS, refSurface);
1625  else {
1626  GlobalTrackerMuonAlignment::trackFitter(itMuon->track(), trackTT, oppositeToMomentum, trackFittedTSOS);
1627  if (trackFittedTSOS.isValid())
1628  extrapolationT = oppositeSmPr.propagate(trackFittedTSOS, refSurface);
1629  }
1630 
1631  if (!extrapolationT.isValid()) {
1632  if (false & alarm)
1633  std::cout << " !!!!!!!!! Catastrophe Out-In extrapolationT.isValid "
1634  << "\a\a\a\a\a\a\a\a" << extrapolationT.isValid() << std::endl;
1635  continue;
1636  }
1637  if (debug_)
1638  std::cout << " extrapolationT.isValid " << extrapolationT.isValid() << std::endl;
1639 
1640  tsosMuonIf = 2;
1641  Rt = GlobalVector((extrapolationT.globalPosition()).x(),
1642  (extrapolationT.globalPosition()).y(),
1643  (extrapolationT.globalPosition()).z());
1644 
1645  Pt = extrapolationT.globalMomentum();
1646  // global parameters of muon
1647  GRm = GlobalVector(
1648  (outerMuTSOS.globalPosition()).x(), (outerMuTSOS.globalPosition()).y(), (outerMuTSOS.globalPosition()).z());
1649  GPm = outerMuTSOS.globalMomentum();
1650  Rt0 = GlobalVector((innerTrackTSOS.globalPosition()).x(),
1651  (innerTrackTSOS.globalPosition()).y(),
1652  (innerTrackTSOS.globalPosition()).z());
1653  Cm = AlgebraicSymMatrix66(extrapolationT.cartesianError().matrix() + outerMuTSOS.cartesianError().matrix());
1654  C0 = AlgebraicSymMatrix66(innerTrackTSOS.cartesianError().matrix());
1655  Ce = AlgebraicSymMatrix66(extrapolationT.cartesianError().matrix());
1656  C1 = AlgebraicSymMatrix66(outerMuTSOS.cartesianError().matrix());
1657 
1658  if (refitMuon_)
1659  GlobalTrackerMuonAlignment::muonFitter(itMuon->outerTrack(), muTT, alongMomentum, muonFittedTSOS);
1660 
1661  if (false & debug_) {
1662  //std::cout<<" ->propDir "<<propagator->propagationDirection()<<std::endl;
1663  PropagationDirectionChooser Chooser;
1664  std::cout << " propDirCh " << Chooser.operator()(*innerTrackTSOS.freeState(), refSurface)
1665  << " Ch == oppisite "
1666  << (oppositeToMomentum == Chooser.operator()(*innerTrackTSOS.freeState(), refSurface)) << std::endl;
1667  std::cout << " --- Nlocal " << Nl.x() << " " << Nl.y() << " " << Nl.z() << " "
1668  << " alfa " << int(asin(Nl.x()) * 57.29578) << std::endl;
1669  }
1670  } // enf of ---- In - Out -----
1671 
1672  else if ((distanceOutOut <= distanceInOut) & (distanceOutOut <= distanceInIn) &
1673  (distanceOutOut <= distanceOutIn)) { // ----- Out - Out ------
1674  if (debug_)
1675  std::cout << " ----- Out - Out -----" << std::endl;
1676 
1677  // reject: momentum of track has opposite direction to muon track
1678  continue;
1679 
1680  const Surface& refSurface = outerMuTSOS.surface();
1681  ConstReferenceCountingPointer<TangentPlane> tpMuGlobal(refSurface.tangentPlane(outerMuTSOS.globalPosition()));
1682  Nl = tpMuGlobal->normalVector();
1683 
1684  // extrapolation to outermost muon hit
1685  extrapolationT = alongSmPr.propagate(outerTrackTSOS, refSurface);
1686 
1687  if (!extrapolationT.isValid()) {
1688  if (alarm)
1689  std::cout << " !!!!!!!!! Catastrophe Out-Out extrapolationT.isValid "
1690  << "\a\a\a\a\a\a\a\a" << extrapolationT.isValid() << std::endl;
1691  continue;
1692  }
1693  if (debug_)
1694  std::cout << " extrapolationT.isValid " << extrapolationT.isValid() << std::endl;
1695 
1696  tsosMuonIf = 2;
1697  Rt = GlobalVector((extrapolationT.globalPosition()).x(),
1698  (extrapolationT.globalPosition()).y(),
1699  (extrapolationT.globalPosition()).z());
1700 
1701  Pt = extrapolationT.globalMomentum();
1702  // global parameters of muon
1703  GRm = GlobalVector(
1704  (outerMuTSOS.globalPosition()).x(), (outerMuTSOS.globalPosition()).y(), (outerMuTSOS.globalPosition()).z());
1705  GPm = outerMuTSOS.globalMomentum();
1706  Rt0 = GlobalVector((outerTrackTSOS.globalPosition()).x(),
1707  (outerTrackTSOS.globalPosition()).y(),
1708  (outerTrackTSOS.globalPosition()).z());
1709  Cm = AlgebraicSymMatrix66(extrapolationT.cartesianError().matrix() + outerMuTSOS.cartesianError().matrix());
1710  C0 = AlgebraicSymMatrix66(outerTrackTSOS.cartesianError().matrix());
1711  Ce = AlgebraicSymMatrix66(extrapolationT.cartesianError().matrix());
1712  C1 = AlgebraicSymMatrix66(outerMuTSOS.cartesianError().matrix());
1713 
1714  if (debug_) {
1715  PropagationDirectionChooser Chooser;
1716  std::cout << " propDirCh " << Chooser.operator()(*outerTrackTSOS.freeState(), refSurface) << " Ch == along "
1717  << (alongMomentum == Chooser.operator()(*outerTrackTSOS.freeState(), refSurface)) << std::endl;
1718  std::cout << " --- Nlocal " << Nl.x() << " " << Nl.y() << " " << Nl.z() << " "
1719  << " alfa " << int(asin(Nl.x()) * 57.29578) << std::endl;
1720  std::cout << " Nornal " << Nl.x() << " " << Nl.y() << " " << Nl.z() << " "
1721  << " alfa " << int(asin(Nl.x()) * 57.29578) << std::endl;
1722  }
1723  } // enf of ---- Out - Out -----
1724  else {
1725  if (alarm)
1726  std::cout << " ------- !!!!!!!!!!!!!!! No proper Out - In \a\a\a\a\a\a\a" << std::endl;
1727  continue;
1728  }
1729 
1730  } // ------------------------------- end Cosmic Muon -----
1731 
1732  if (tsosMuonIf == 0) {
1733  if (info) {
1734  std::cout << "No tsosMuon !!!!!!" << std::endl;
1735  continue;
1736  }
1737  }
1738  TrajectoryStateOnSurface tsosMuon;
1739  if (tsosMuonIf == 1)
1740  tsosMuon = muTT.innermostMeasurementState();
1741  else
1742  tsosMuon = muTT.outermostMeasurementState();
1743 
1744  //GlobalTrackerMuonAlignment::misalignMuon(GRm, GPm, Nl, Rt, Rm, Pm);
1745  AlgebraicVector4 LPRm; // muon local (dx/dz, dy/dz, x, y)
1746  GlobalTrackerMuonAlignment::misalignMuonL(GRm, GPm, Nl, Rt, Rm, Pm, LPRm, extrapolationT, tsosMuon);
1747 
1748  if (refitTrack_) {
1749  if (!trackFittedTSOS.isValid()) {
1750  if (info)
1751  std::cout << " ================= trackFittedTSOS notValid !!!!!!!! " << std::endl;
1752  continue;
1753  }
1754  if (debug_)
1755  this->debugTrajectorySOS(" trackFittedTSOS ", trackFittedTSOS);
1756  }
1757 
1758  if (refitMuon_) {
1759  if (!muonFittedTSOS.isValid()) {
1760  if (info)
1761  std::cout << " ================= muonFittedTSOS notValid !!!!!!!! " << std::endl;
1762  continue;
1763  }
1764  if (debug_)
1765  this->debugTrajectorySOS(" muonFittedTSOS ", muonFittedTSOS);
1766  Rm = GlobalVector((muonFittedTSOS.globalPosition()).x(),
1767  (muonFittedTSOS.globalPosition()).y(),
1768  (muonFittedTSOS.globalPosition()).z());
1769  Pm = muonFittedTSOS.globalMomentum();
1770  LPRm = AlgebraicVector4(muonFittedTSOS.localParameters().vector()(1),
1771  muonFittedTSOS.localParameters().vector()(2),
1772  muonFittedTSOS.localParameters().vector()(3),
1773  muonFittedTSOS.localParameters().vector()(4));
1774  }
1775  GlobalVector resR = Rm - Rt;
1776  GlobalVector resP0 = Pm - Pt;
1777  GlobalVector resP = Pm / Pm.mag() - Pt / Pt.mag();
1778  float RelMomResidual = (Pm.mag() - Pt.mag()) / (Pt.mag() + 1.e-6);
1779  ;
1780 
1781  AlgebraicVector6 Vm;
1782  Vm(0) = resR.x();
1783  Vm(1) = resR.y();
1784  Vm(2) = resR.z();
1785  Vm(3) = resP0.x();
1786  Vm(4) = resP0.y();
1787  Vm(5) = resP0.z();
1788  float Rmuon = Rm.perp();
1789  float Zmuon = Rm.z();
1790  float alfa_x = atan2(Nl.x(), Nl.y()) * 57.29578;
1791 
1792  if (debug_) {
1793  std::cout << " Nx Ny Nz alfa_x " << Nl.x() << " " << Nl.y() << " " << Nl.z() << " " << alfa_x << std::endl;
1794  //std::cout<<" Rm "<<Rm<<std::endl<<" Rt "<<Rt<<std::endl<<" resR "<<resR<<std::endl
1795  // <<" resP "<<resP<<" dp/p "<<RelMomResidual<<std::endl;
1796  }
1797 
1798  double chi_d = 0;
1799  for (int i = 0; i <= 5; i++)
1800  chi_d += Vm(i) * Vm(i) / Cm(i, i);
1801 
1802  AlgebraicVector5 Vml(tsosMuon.localParameters().vector() - extrapolationT.localParameters().vector());
1803  AlgebraicSymMatrix55 m(tsosMuon.localError().matrix() + extrapolationT.localError().matrix());
1804  AlgebraicSymMatrix55 Cml(tsosMuon.localError().matrix() + extrapolationT.localError().matrix());
1805  bool ierrLoc = !m.Invert();
1806  if (ierrLoc && debug_ && info) {
1807  std::cout << " ==== Error inverting Local covariance matrix ==== " << std::endl;
1808  continue;
1809  }
1810  double chi_Loc = ROOT::Math::Similarity(Vml, m);
1811  if (debug_)
1812  std::cout << " chi_Loc px/pz/err " << chi_Loc << " " << Vml(1) / sqrt(Cml(1, 1)) << std::endl;
1813 
1814  if (Pt.mag() < 15.)
1815  continue;
1816  if (Pm.mag() < 5.)
1817  continue;
1818 
1819  //if(Pt.mag() < 30.) continue; // momenum cut < 30GeV
1820  //if(Pt.mag() < 60.) continue; // momenum cut < 30GeV
1821  //if(Pt.mag() > 50.) continue; // momenum cut > 50GeV
1822  //if(Pt.mag() > 100.) continue; // momenum cut > 100GeV
1823  //if(trackTT.charge() < 0) continue; // select positive charge
1824  //if(trackTT.charge() > 0) continue; // select negative charge
1825 
1826  //if(fabs(resR.x()) > 5.) continue; // strong cut X
1827  //if(fabs(resR.y()) > 5.) continue; // Y
1828  //if(fabs(resR.z()) > 5.) continue; // Z
1829  //if(fabs(resR.mag()) > 7.5) continue; // dR
1830 
1831  //if(fabs(RelMomResidual) > 0.5) continue;
1832  if (fabs(resR.x()) > 20.)
1833  continue;
1834  if (fabs(resR.y()) > 20.)
1835  continue;
1836  if (fabs(resR.z()) > 20.)
1837  continue;
1838  if (fabs(resR.mag()) > 30.)
1839  continue;
1840  if (fabs(resP.x()) > 0.06)
1841  continue;
1842  if (fabs(resP.y()) > 0.06)
1843  continue;
1844  if (fabs(resP.z()) > 0.06)
1845  continue;
1846  if (chi_d > 40.)
1847  continue;
1848 
1849  // select Barrel
1850  //if(Rmuon < 400. || Rmuon > 450.) continue;
1851  //if(Zmuon < -600. || Zmuon > 600.) continue;
1852  //if(fabs(Nl.z()) > 0.95) continue;
1853  //MuSelect = " Barrel";
1854  // EndCap1
1855  //if(Rmuon < 120. || Rmuon > 450.) continue;
1856  //if(Zmuon < -720.) continue;
1857  //if(Zmuon > -580.) continue;
1858  //if(fabs(Nl.z()) < 0.95) continue;
1859  //MuSelect = " EndCap1";
1860  // EndCap2
1861  //if(Rmuon < 120. || Rmuon > 450.) continue;
1862  //if(Zmuon > 720.) continue;
1863  //if(Zmuon < 580.) continue;
1864  //if(fabs(Nl.z()) < 0.95) continue;
1865  //MuSelect = " EndCap2";
1866  // select All
1867  if (Rmuon < 120. || Rmuon > 450.)
1868  continue;
1869  if (Zmuon < -720. || Zmuon > 720.)
1870  continue;
1871  MuSelect = " Barrel+EndCaps";
1872 
1873  if (debug_)
1874  std::cout << " .............. passed all cuts" << std::endl;
1875 
1876  N_track++;
1877  // gradient and Hessian for each track
1878 
1879  GlobalTrackerMuonAlignment::gradientGlobalAlg(Rt, Pt, Rm, Nl, Cm);
1880  GlobalTrackerMuonAlignment::gradientGlobal(Rt, Pt, Rm, Pm, Nl, Cm);
1881 
1882  CLHEP::HepSymMatrix covLoc(4, 0);
1883  for (int i = 1; i <= 4; i++)
1884  for (int j = 1; j <= i; j++) {
1885  covLoc(i, j) = (tsosMuon.localError().matrix() + extrapolationT.localError().matrix())(i, j);
1886  //if(i != j) Cov(i,j) = 0.;
1887  }
1888 
1889  const Surface& refSurface = tsosMuon.surface();
1890  CLHEP::HepMatrix rotLoc(3, 3, 0);
1891  rotLoc(1, 1) = refSurface.rotation().xx();
1892  rotLoc(1, 2) = refSurface.rotation().xy();
1893  rotLoc(1, 3) = refSurface.rotation().xz();
1894 
1895  rotLoc(2, 1) = refSurface.rotation().yx();
1896  rotLoc(2, 2) = refSurface.rotation().yy();
1897  rotLoc(2, 3) = refSurface.rotation().yz();
1898 
1899  rotLoc(3, 1) = refSurface.rotation().zx();
1900  rotLoc(3, 2) = refSurface.rotation().zy();
1901  rotLoc(3, 3) = refSurface.rotation().zz();
1902 
1903  CLHEP::HepVector posLoc(3);
1904  posLoc(1) = refSurface.position().x();
1905  posLoc(2) = refSurface.position().y();
1906  posLoc(3) = refSurface.position().z();
1907 
1908  GlobalTrackerMuonAlignment::gradientLocal(Rt, Pt, Rm, Pm, Nl, covLoc, rotLoc, posLoc, LPRm);
1909 
1910  if (debug_) {
1911  std::cout << " Norm " << Nl << std::endl;
1912  std::cout << " posLoc " << posLoc.T() << std::endl;
1913  std::cout << " rotLoc " << rotLoc << std::endl;
1914  }
1915 
1916  // ----------------------------------------------------- fill histogram
1917  histo->Fill(itMuon->track()->pt());
1918 
1919  //histo2->Fill(itMuon->track()->outerP());
1920  histo2->Fill(Pt.mag());
1921  histo3->Fill((PI / 2. - itMuon->track()->outerTheta()));
1922  histo4->Fill(itMuon->track()->phi());
1923  histo5->Fill(Rmuon);
1924  histo6->Fill(Zmuon);
1925  histo7->Fill(RelMomResidual);
1926  //histo8->Fill(chi);
1927  histo8->Fill(chi_d);
1928 
1929  histo101->Fill(Zmuon, Rmuon);
1930  histo101->Fill(Rt0.z(), Rt0.perp());
1931  histo102->Fill(Rt0.x(), Rt0.y());
1932  histo102->Fill(Rm.x(), Rm.y());
1933 
1934  histo11->Fill(resR.mag());
1935  if (fabs(Nl.x()) < 0.98)
1936  histo12->Fill(resR.x());
1937  if (fabs(Nl.y()) < 0.98)
1938  histo13->Fill(resR.y());
1939  if (fabs(Nl.z()) < 0.98)
1940  histo14->Fill(resR.z());
1941  histo15->Fill(resP.x());
1942  histo16->Fill(resP.y());
1943  histo17->Fill(resP.z());
1944 
1945  if ((fabs(Nl.x()) < 0.98) && (fabs(Nl.y()) < 0.98) && (fabs(Nl.z()) < 0.98)) {
1946  histo18->Fill(sqrt(C0(0, 0)));
1947  histo19->Fill(sqrt(C1(0, 0)));
1948  histo20->Fill(sqrt(C1(0, 0) + Ce(0, 0)));
1949  }
1950  if (fabs(Nl.x()) < 0.98)
1951  histo21->Fill(Vm(0) / sqrt(Cm(0, 0)));
1952  if (fabs(Nl.y()) < 0.98)
1953  histo22->Fill(Vm(1) / sqrt(Cm(1, 1)));
1954  if (fabs(Nl.z()) < 0.98)
1955  histo23->Fill(Vm(2) / sqrt(Cm(2, 2)));
1956  histo24->Fill(Vm(3) / sqrt(C1(3, 3) + Ce(3, 3)));
1957  histo25->Fill(Vm(4) / sqrt(C1(4, 4) + Ce(4, 4)));
1958  histo26->Fill(Vm(5) / sqrt(C1(5, 5) + Ce(5, 5)));
1959  histo27->Fill(Nl.x());
1960  histo28->Fill(Nl.y());
1961  histo29->Fill(lenghtTrack);
1962  histo30->Fill(lenghtMuon);
1963  histo31->Fill(chi_Loc);
1964  histo32->Fill(Vml(1) / sqrt(Cml(1, 1)));
1965  histo33->Fill(Vml(2) / sqrt(Cml(2, 2)));
1966  histo34->Fill(Vml(3) / sqrt(Cml(3, 3)));
1967  histo35->Fill(Vml(4) / sqrt(Cml(4, 4)));
1968 
1969  if (debug_) { //--------------------------------- debug print ----------
1970 
1971  std::cout << " diag 0[ " << C0(0, 0) << " " << C0(1, 1) << " " << C0(2, 2) << " " << C0(3, 3) << " " << C0(4, 4)
1972  << " " << C0(5, 5) << " ]" << std::endl;
1973  std::cout << " diag e[ " << Ce(0, 0) << " " << Ce(1, 1) << " " << Ce(2, 2) << " " << Ce(3, 3) << " " << Ce(4, 4)
1974  << " " << Ce(5, 5) << " ]" << std::endl;
1975  std::cout << " diag 1[ " << C1(0, 0) << " " << C1(1, 1) << " " << C1(2, 2) << " " << C1(3, 3) << " " << C1(4, 4)
1976  << " " << C1(5, 5) << " ]" << std::endl;
1977  std::cout << " Rm " << Rm.x() << " " << Rm.y() << " " << Rm.z() << " Pm " << Pm.x() << " " << Pm.y() << " "
1978  << Pm.z() << std::endl;
1979  std::cout << " Rt " << Rt.x() << " " << Rt.y() << " " << Rt.z() << " Pt " << Pt.x() << " " << Pt.y() << " "
1980  << Pt.z() << std::endl;
1981  std::cout << " Nl*(Rm-Rt) " << Nl.dot(Rm - Rt) << std::endl;
1982  std::cout << " resR " << resR.x() << " " << resR.y() << " " << resR.z() << " resP " << resP.x() << " " << resP.y()
1983  << " " << resP.z() << std::endl;
1984  std::cout << " Rm-t " << (Rm - Rt).x() << " " << (Rm - Rt).y() << " " << (Rm - Rt).z() << " Pm-t "
1985  << (Pm - Pt).x() << " " << (Pm - Pt).y() << " " << (Pm - Pt).z() << std::endl;
1986  std::cout << " Vm " << Vm << std::endl;
1987  std::cout << " +- " << sqrt(Cm(0, 0)) << " " << sqrt(Cm(1, 1)) << " " << sqrt(Cm(2, 2)) << " " << sqrt(Cm(3, 3))
1988  << " " << sqrt(Cm(4, 4)) << " " << sqrt(Cm(5, 5)) << std::endl;
1989  std::cout << " Pmuon Ptrack dP/Ptrack " << itMuon->outerTrack()->p() << " " << itMuon->track()->outerP() << " "
1990  << (itMuon->outerTrack()->p() - itMuon->track()->outerP()) / itMuon->track()->outerP() << std::endl;
1991  std::cout << " cov matrix " << std::endl;
1992  std::cout << Cm << std::endl;
1993  std::cout << " diag [ " << Cm(0, 0) << " " << Cm(1, 1) << " " << Cm(2, 2) << " " << Cm(3, 3) << " " << Cm(4, 4)
1994  << " " << Cm(5, 5) << " ]" << std::endl;
1995 
1996  AlgebraicSymMatrix66 Ro;
1997  double Diag[6];
1998  for (int i = 0; i <= 5; i++)
1999  Diag[i] = sqrt(Cm(i, i));
2000  for (int i = 0; i <= 5; i++)
2001  for (int j = 0; j <= 5; j++)
2002  Ro(i, j) = Cm(i, j) / Diag[i] / Diag[j];
2003  std::cout << " correlation matrix " << std::endl;
2004  std::cout << Ro << std::endl;
2005 
2006  AlgebraicSymMatrix66 CmI;
2007  for (int i = 0; i <= 5; i++)
2008  for (int j = 0; j <= 5; j++)
2009  CmI(i, j) = Cm(i, j);
2010 
2011  bool ierr = !CmI.Invert();
2012  if (ierr) {
2013  if (alarm || debug_)
2014  std::cout << " Error inverse covariance matrix !!!!!!!!!!!" << std::endl;
2015  continue;
2016  }
2017  std::cout << " inverse cov matrix " << std::endl;
2018  std::cout << Cm << std::endl;
2019 
2020  double chi = ROOT::Math::Similarity(Vm, CmI);
2021  std::cout << " chi chi_d " << chi << " " << chi_d << std::endl;
2022  } // end of debug_ printout --------------------------------------------
2023 
2024  } // end loop on selected muons, i.e. Jim's globalMuon
2025 
2026 } //end of analyzeTrackTrajectory
GlobalTrackerMuonAlignment::trackingGeometry_
const GlobalTrackingGeometry * trackingGeometry_
Definition: GlobalTrackerMuonAlignment.cc:182
HLT_2022v11_cff.magneticField
magneticField
Definition: HLT_2022v11_cff.py:487
GlobalTrackerMuonAlignment::watchGlobalPositionRcd_
edm::ESWatcher< GlobalPositionRcd > watchGlobalPositionRcd_
Definition: GlobalTrackerMuonAlignment.cc:189
GlobalTrackerMuonAlignment::theTrackingGeometry
edm::ESHandle< GlobalTrackingGeometry > theTrackingGeometry
Definition: GlobalTrackerMuonAlignment.cc:181
info
static const TGPicture * info(bool iBackgroundIsBlack)
Definition: FWCollectionSummaryWidget.cc:153
PV3DBase::perp
T perp() const
Definition: PV3DBase.h:69
mps_fire.i
i
Definition: mps_fire.py:428
SmartPropagator::propagate
TrajectoryStateOnSurface propagate(STA const &state, SUR const &surface) const
Definition: Propagator.h:50
GlobalTrackerMuonAlignment::gradientGlobalAlg
void gradientGlobalAlg(GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, AlgebraicSymMatrix66 &)
Definition: GlobalTrackerMuonAlignment.cc:2029
GlobalTrackerMuonAlignment::gradientLocal
void gradientLocal(GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, CLHEP::HepSymMatrix &, CLHEP::HepMatrix &, CLHEP::HepVector &, AlgebraicVector4 &)
Definition: GlobalTrackerMuonAlignment.cc:2305
TkRotation::xx
T xx() const
Definition: extTkRotation.h:193
TrajectoryStateOnSurface::localError
const LocalTrajectoryError & localError() const
Definition: TrajectoryStateOnSurface.h:77
Surface::tangentPlane
virtual ConstReferenceCountingPointer< TangentPlane > tangentPlane(const GlobalPoint &) const =0
PV3DBase::z
T z() const
Definition: PV3DBase.h:61
TrackCandidateProducer_cfi.propagator
propagator
Definition: TrackCandidateProducer_cfi.py:17
Vector3DBase::dot
PreciseFloatType< T, U >::Type dot(const Vector3DBase< U, FrameTag > &v) const
Definition: Vector3DBase.h:99
TkRotation::xy
T xy() const
Definition: extTkRotation.h:194
GlobalTrackerMuonAlignment::m_globalPosToken
edm::ESGetToken< Alignments, GlobalPositionRcd > m_globalPosToken
Definition: GlobalTrackerMuonAlignment.cc:160
GlobalTrackerMuonAlignment::m_TkGeometryToken
edm::ESGetToken< GlobalTrackingGeometry, GlobalTrackingGeometryRecord > m_TkGeometryToken
Definition: GlobalTrackerMuonAlignment.cc:158
SteppingHelixPropagator_cfi.SteppingHelixPropagator
SteppingHelixPropagator
Definition: SteppingHelixPropagator_cfi.py:3
TkRotation::zz
T zz() const
Definition: extTkRotation.h:201
GlobalTrackerMuonAlignment::iteratorHcalRcd
std::vector< AlignTransform >::const_iterator iteratorHcalRcd
Definition: GlobalTrackerMuonAlignment.cc:217
TkRotation::yy
T yy() const
Definition: extTkRotation.h:197
TrajectoryStateOnSurface::localParameters
const LocalTrajectoryParameters & localParameters() const
Definition: TrajectoryStateOnSurface.h:73
TkRotation::yz
T yz() const
Definition: extTkRotation.h:198
MillePedeFileConverter_cfg.e
e
Definition: MillePedeFileConverter_cfg.py:37
TrajectoryStateOnSurface::surface
const SurfaceType & surface() const
Definition: TrajectoryStateOnSurface.h:78
GlobalTrackerMuonAlignment::gradientGlobal
void gradientGlobal(GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, AlgebraicSymMatrix66 &)
Definition: GlobalTrackerMuonAlignment.cc:2120
Alignments::m_align
std::vector< AlignTransform > m_align
Definition: Alignments.h:19
Plane
Definition: Plane.h:16
TrajectoryStateOnSurface::cartesianError
const CartesianTrajectoryError cartesianError() const
Definition: TrajectoryStateOnSurface.h:71
TkRotation::zx
T zx() const
Definition: extTkRotation.h:199
PV3DBase::x
T x() const
Definition: PV3DBase.h:59
PV3DBase::y
T y() const
Definition: PV3DBase.h:60
Chi2MeasurementEstimator_cfi.Chi2MeasurementEstimator
Chi2MeasurementEstimator
Definition: Chi2MeasurementEstimator_cfi.py:5
iEvent
int iEvent
Definition: GenABIO.cc:224
TrajectoryStateOnSurface::globalPosition
GlobalPoint globalPosition() const
Definition: TrajectoryStateOnSurface.h:65
GlobalTrackerMuonAlignment::watchTrackingGeometry_
edm::ESWatcher< GlobalTrackingGeometryRecord > watchTrackingGeometry_
Definition: GlobalTrackerMuonAlignment.cc:180
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
LocalTrajectoryParameters::vector
AlgebraicVector5 vector() const
Definition: LocalTrajectoryParameters.h:120
visualization-live-secondInstance_cfg.m
m
Definition: visualization-live-secondInstance_cfg.py:80
TkRotation::zy
T zy() const
Definition: extTkRotation.h:200
PV3DBase::mag
T mag() const
Definition: PV3DBase.h:64
GlobalTrackerMuonAlignment::iteratorMuonRcd
std::vector< AlignTransform >::const_iterator iteratorMuonRcd
Definition: GlobalTrackerMuonAlignment.cc:215
AlgebraicVector5
ROOT::Math::SVector< double, 5 > AlgebraicVector5
Definition: AlgebraicROOTObjects.h:14
AlgebraicVector4
ROOT::Math::SVector< double, 4 > AlgebraicVector4
Definition: AlgebraicROOTObjects.h:13
edm::EventSetup::getHandle
ESHandle< T > getHandle(const ESGetToken< T, R > &iToken) const
Definition: EventSetup.h:151
GlobalTrackerMuonAlignment::m_ttrhBuilderToken
edm::ESGetToken< TransientTrackingRecHitBuilder, TransientRecHitRecord > m_ttrhBuilderToken
Definition: GlobalTrackerMuonAlignment.cc:162
edm::ESHandleBase::isValid
bool isValid() const
Definition: ESHandle.h:44
GlobalTrackerMuonAlignment::misalignMuonL
void misalignMuonL(GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, GlobalVector &, AlgebraicVector4 &, TrajectoryStateOnSurface &, TrajectoryStateOnSurface &)
Definition: GlobalTrackerMuonAlignment.cc:2712
AlgebraicSymMatrix66
ROOT::Math::SMatrix< double, 6, 6, ROOT::Math::MatRepSym< double, 6 > > AlgebraicSymMatrix66
Definition: AlgebraicROOTObjects.h:24
GlobalTrackerMuonAlignment::iteratorEcalRcd
std::vector< AlignTransform >::const_iterator iteratorEcalRcd
Definition: GlobalTrackerMuonAlignment.cc:216
DetId
Definition: DetId.h:17
GlobalTrackerMuonAlignment::trackFitter
void trackFitter(reco::TrackRef, reco::TransientTrack &, PropagationDirection, TrajectoryStateOnSurface &)
Definition: GlobalTrackerMuonAlignment.cc:2994
mag
T mag() const
The vector magnitude. Equivalent to sqrt(vec.mag2())
Definition: Basic3DVectorLD.h:127
tracks
auto const & tracks
cannot be loose
Definition: CAHitNtupletGeneratorKernelsImpl.h:148
GloballyPositioned::position
const PositionType & position() const
Definition: GloballyPositioned.h:36
DetId::rawId
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:57
GlobalTrackerMuonAlignment::TTRHBuilder
const TransientTrackingRecHitBuilder * TTRHBuilder
Definition: GlobalTrackerMuonAlignment.cc:198
SmartPropagator_cfi.SmartPropagator
SmartPropagator
Definition: SmartPropagator_cfi.py:3
AlgebraicSymMatrix55
ROOT::Math::SMatrix< double, 5, 5, ROOT::Math::MatRepSym< double, 5 > > AlgebraicSymMatrix55
Definition: AlgebraicROOTObjects.h:23
GlobalTrackerMuonAlignment::watchMagneticFieldRecord_
edm::ESWatcher< IdealMagneticFieldRecord > watchMagneticFieldRecord_
Definition: GlobalTrackerMuonAlignment.cc:184
TrajectoryStateOnSurface::globalMomentum
GlobalVector globalMomentum() const
Definition: TrajectoryStateOnSurface.h:66
edm::ESWatcher::check
bool check(const edm::EventSetup &iSetup)
Definition: ESWatcher.h:57
GlobalTrackerMuonAlignment::m_MagFieldToken
edm::ESGetToken< MagneticField, IdealMagneticFieldRecord > m_MagFieldToken
Definition: GlobalTrackerMuonAlignment.cc:159
reco
fixed size matrix
Definition: AlignmentAlgorithmBase.h:46
edm
HLT enums.
Definition: AlignableModifier.h:19
LocalTrajectoryError::matrix
const AlgebraicSymMatrix55 & matrix() const
Definition: LocalTrajectoryError.h:60
CartesianTrajectoryError::matrix
const AlgebraicSymMatrix66 & matrix() const
Definition: CartesianTrajectoryError.h:28
TrajectoryStateOnSurface::freeState
FreeTrajectoryState const * freeState(bool withErrors=true) const
Definition: TrajectoryStateOnSurface.h:58
AlgebraicVector6
ROOT::Math::SVector< double, 6 > AlgebraicVector6
Definition: AlgebraicROOTObjects.h:15
Plane::normalVector
GlobalVector normalVector() const
Definition: Plane.h:41
GlobalTrackerMuonAlignment::iteratorTrackerRcd
std::vector< AlignTransform >::const_iterator iteratorTrackerRcd
Definition: GlobalTrackerMuonAlignment.cc:214
GloballyPositioned::rotation
const RotationType & rotation() const
Definition: GloballyPositioned.h:38
PDWG_BPHSkim_cff.muons
muons
Definition: PDWG_BPHSkim_cff.py:47
gather_cfg.cout
cout
Definition: gather_cfg.py:144
GlobalTrackerMuonAlignment::m_propToken
edm::ESGetToken< Propagator, TrackingComponentsRecord > m_propToken
Definition: GlobalTrackerMuonAlignment.cc:161
TkRotation::xz
T xz() const
Definition: extTkRotation.h:195
GlobalTrackerMuonAlignment::muonFitter
void muonFitter(reco::TrackRef, reco::TransientTrack &, PropagationDirection, TrajectoryStateOnSurface &)
Definition: GlobalTrackerMuonAlignment.cc:3109
TkRotation::yx
T yx() const
Definition: extTkRotation.h:196
GlobalTrackerMuonAlignment::MuRHBuilder
MuonTransientTrackingRecHitBuilder * MuRHBuilder
Definition: GlobalTrackerMuonAlignment.cc:197
GlobalVector
Global3DVector GlobalVector
Definition: GlobalVector.h:10
GlobalTrackerMuonAlignment::debugTrajectorySOS
void debugTrajectorySOS(const std::string, TrajectoryStateOnSurface &)
Definition: GlobalTrackerMuonAlignment.cc:3247

◆ beginJob()

void GlobalTrackerMuonAlignment::beginJob ( )
overrideprivatevirtual

Reimplemented from edm::one::EDAnalyzerBase.

Definition at line 324 of file GlobalTrackerMuonAlignment.cc.

References bookHist(), collectionCosmic, collectionIsolated, cosmicMuonMode_, file, Gfr, Grad, GradL, Hess, HessL, mps_fire::i, Inf, isolatedMuonMode_, dqmiolumiharvest::j, N_event, N_track, and rootOutFile_.

324  {
325  N_event = 0;
326  N_track = 0;
327 
328  if (cosmicMuonMode_ == true && isolatedMuonMode_ == false) {
329  collectionCosmic = true;
330  collectionIsolated = false;
331  } else if (cosmicMuonMode_ == false && isolatedMuonMode_ == true) {
332  collectionCosmic = false;
333  collectionIsolated = true;
334  } else {
335  collectionCosmic = false;
336  collectionIsolated = false;
337  }
338 
339  for (int i = 0; i <= 2; i++) {
340  Gfr(i) = 0.;
341  for (int j = 0; j <= 2; j++) {
342  Inf(i, j) = 0.;
343  }
344  }
345 
346  Grad = CLHEP::HepVector(6, 0);
347  Hess = CLHEP::HepSymMatrix(6, 0);
348 
349  GradL = CLHEP::HepVector(6, 0);
350  HessL = CLHEP::HepSymMatrix(6, 0);
351 
352  // histograms
353  TDirectory* dirsave = gDirectory;
354 
355  file = new TFile(rootOutFile_.c_str(), "recreate");
356  const bool oldAddDir = TH1::AddDirectoryStatus();
357 
358  TH1::AddDirectory(true);
359 
360  this->bookHist();
361 
362  TH1::AddDirectory(oldAddDir);
363  dirsave->cd();
364 }
mps_fire.i
i
Definition: mps_fire.py:428

◆ bookHist()

void GlobalTrackerMuonAlignment::bookHist ( )

Definition at line 477 of file GlobalTrackerMuonAlignment.cc.

References histo, histo101, histo102, histo11, histo12, histo13, histo14, histo15, histo16, histo17, histo18, histo19, histo2, histo20, histo21, histo22, histo23, histo24, histo25, histo26, histo27, histo28, histo29, histo3, histo30, histo31, histo32, histo33, histo34, histo35, histo4, histo5, histo6, histo7, and histo8.

Referenced by beginJob().

477  {
478  double PI = 3.1415927;
479  histo = new TH1F("Pt", "pt", 1000, 0, 100);
480  histo2 = new TH1F("P", "P [GeV/c]", 400, 0., 400.);
481  histo2->GetXaxis()->SetTitle("momentum [GeV/c]");
482  histo3 = new TH1F("outerLambda", "#lambda outer", 100, -PI / 2., PI / 2.);
483  histo3->GetXaxis()->SetTitle("#lambda outer");
484  histo4 = new TH1F("phi", "#phi [rad]", 100, -PI, PI);
485  histo4->GetXaxis()->SetTitle("#phi [rad]");
486  histo5 = new TH1F("Rmuon", "inner muon hit R [cm]", 100, 0., 800.);
487  histo5->GetXaxis()->SetTitle("R of muon [cm]");
488  histo6 = new TH1F("Zmuon", "inner muon hit Z[cm]", 100, -1000., 1000.);
489  histo6->GetXaxis()->SetTitle("Z of muon [cm]");
490  histo7 = new TH1F("(Pm-Pt)/Pt", " (Pmuon-Ptrack)/Ptrack", 100, -2., 2.);
491  histo7->GetXaxis()->SetTitle("(Pmuon-Ptrack)/Ptrack");
492  histo8 = new TH1F("chi muon-track", "#chi^{2}(muon-track)", 1000, 0., 1000.);
493  histo8->GetXaxis()->SetTitle("#chi^{2} of muon w.r.t. propagated track");
494  histo11 = new TH1F("distance muon-track", "distance muon w.r.t track [cm]", 100, 0., 30.);
495  histo11->GetXaxis()->SetTitle("distance of muon w.r.t. track [cm]");
496  histo12 = new TH1F("Xmuon-Xtrack", "Xmuon-Xtrack [cm]", 200, -20., 20.);
497  histo12->GetXaxis()->SetTitle("Xmuon - Xtrack [cm]");
498  histo13 = new TH1F("Ymuon-Ytrack", "Ymuon-Ytrack [cm]", 200, -20., 20.);
499  histo13->GetXaxis()->SetTitle("Ymuon - Ytrack [cm]");
500  histo14 = new TH1F("Zmuon-Ztrack", "Zmuon-Ztrack [cm]", 200, -20., 20.);
501  histo14->GetXaxis()->SetTitle("Zmuon-Ztrack [cm]");
502  histo15 = new TH1F("NXmuon-NXtrack", "NXmuon-NXtrack [rad]", 200, -.1, .1);
503  histo15->GetXaxis()->SetTitle("N_{X}(muon)-N_{X}(track) [rad]");
504  histo16 = new TH1F("NYmuon-NYtrack", "NYmuon-NYtrack [rad]", 200, -.1, .1);
505  histo16->GetXaxis()->SetTitle("N_{Y}(muon)-N_{Y}(track) [rad]");
506  histo17 = new TH1F("NZmuon-NZtrack", "NZmuon-NZtrack [rad]", 200, -.1, .1);
507  histo17->GetXaxis()->SetTitle("N_{Z}(muon)-N_{Z}(track) [rad]");
508  histo18 = new TH1F("expected error of Xinner", "outer hit of inner tracker", 100, 0, .01);
509  histo18->GetXaxis()->SetTitle("expected error of Xinner [cm]");
510  histo19 = new TH1F("expected error of Xmuon", "inner hit of muon", 100, 0, .1);
511  histo19->GetXaxis()->SetTitle("expected error of Xmuon [cm]");
512  histo20 = new TH1F("expected error of Xmuon-Xtrack", "muon w.r.t. propagated track", 100, 0., 10.);
513  histo20->GetXaxis()->SetTitle("expected error of Xmuon-Xtrack [cm]");
514  histo21 = new TH1F("pull of Xmuon-Xtrack", "pull of Xmuon-Xtrack", 100, -10., 10.);
515  histo21->GetXaxis()->SetTitle("(Xmuon-Xtrack)/expected error");
516  histo22 = new TH1F("pull of Ymuon-Ytrack", "pull of Ymuon-Ytrack", 100, -10., 10.);
517  histo22->GetXaxis()->SetTitle("(Ymuon-Ytrack)/expected error");
518  histo23 = new TH1F("pull of Zmuon-Ztrack", "pull of Zmuon-Ztrack", 100, -10., 10.);
519  histo23->GetXaxis()->SetTitle("(Zmuon-Ztrack)/expected error");
520  histo24 = new TH1F("pull of PXmuon-PXtrack", "pull of PXmuon-PXtrack", 100, -10., 10.);
521  histo24->GetXaxis()->SetTitle("(P_{X}(muon)-P_{X}(track))/expected error");
522  histo25 = new TH1F("pull of PYmuon-PYtrack", "pull of PYmuon-PYtrack", 100, -10., 10.);
523  histo25->GetXaxis()->SetTitle("(P_{Y}(muon)-P_{Y}(track))/expected error");
524  histo26 = new TH1F("pull of PZmuon-PZtrack", "pull of PZmuon-PZtrack", 100, -10., 10.);
525  histo26->GetXaxis()->SetTitle("(P_{Z}(muon)-P_{Z}(track))/expected error");
526  histo27 = new TH1F("N_x", "Nx of tangent plane", 120, -1.1, 1.1);
527  histo27->GetXaxis()->SetTitle("normal vector projection N_{X}");
528  histo28 = new TH1F("N_y", "Ny of tangent plane", 120, -1.1, 1.1);
529  histo28->GetXaxis()->SetTitle("normal vector projection N_{Y}");
530  histo29 = new TH1F("lenght of track", "lenght of track", 200, 0., 400);
531  histo29->GetXaxis()->SetTitle("lenght of track [cm]");
532  histo30 = new TH1F("lenght of muon", "lenght of muon", 200, 0., 800);
533  histo30->GetXaxis()->SetTitle("lenght of muon [cm]");
534 
535  histo31 = new TH1F("local chi muon-track", "#local chi^{2}(muon-track)", 1000, 0., 1000.);
536  histo31->GetXaxis()->SetTitle("#local chi^{2} of muon w.r.t. propagated track");
537  histo32 = new TH1F("pull of Px/Pz local", "pull of Px/Pz local", 100, -10., 10.);
538  histo32->GetXaxis()->SetTitle("local (Px/Pz(muon) - Px/Pz(track))/expected error");
539  histo33 = new TH1F("pull of Py/Pz local", "pull of Py/Pz local", 100, -10., 10.);
540  histo33->GetXaxis()->SetTitle("local (Py/Pz(muon) - Py/Pz(track))/expected error");
541  histo34 = new TH1F("pull of X local", "pull of X local", 100, -10., 10.);
542  histo34->GetXaxis()->SetTitle("local (Xmuon - Xtrack)/expected error");
543  histo35 = new TH1F("pull of Y local", "pull of Y local", 100, -10., 10.);
544  histo35->GetXaxis()->SetTitle("local (Ymuon - Ytrack)/expected error");
545 
546  histo101 = new TH2F("Rtr/mu vs Ztr/mu", "hit of track/muon", 100, -800., 800., 100, 0., 600.);
547  histo101->GetXaxis()->SetTitle("Z of track/muon [cm]");
548  histo101->GetYaxis()->SetTitle("R of track/muon [cm]");
549  histo102 = new TH2F("Ytr/mu vs Xtr/mu", "hit of track/muon", 100, -600., 600., 100, -600., 600.);
550  histo102->GetXaxis()->SetTitle("X of track/muon [cm]");
551  histo102->GetYaxis()->SetTitle("Y of track/muon [cm]");
552 }

◆ CLHEP_dot()

double GlobalTrackerMuonAlignment::CLHEP_dot ( const CLHEP::HepVector &  a,
const CLHEP::HepVector &  b 
)
inline

Definition at line 148 of file GlobalTrackerMuonAlignment.cc.

References a, and b.

Referenced by gradientGlobal(), gradientLocal(), misalignMuon(), and misalignMuonL().

148  {
149  return a(1) * b(1) + a(2) * b(2) + a(3) * b(3);
150  }
b
double b
Definition: hdecay.h:118
a
double a
Definition: hdecay.h:119

◆ debugTrackHit() [1/2]

void GlobalTrackerMuonAlignment::debugTrackHit ( const std::string  title,
reco::TrackRef  alongTr 
)

Definition at line 3228 of file GlobalTrackerMuonAlignment.cc.

References gather_cfg::cout, DetId::Muon, runGCPTkAlMap::title, and DetId::Tracker.

Referenced by muonFitter(), and trackFitter().

3228  {
3229  std::cout << " ------- " << title << " --------" << std::endl;
3230  int nHit = 1;
3231  for (auto const& hit : alongTr->recHits()) {
3232  std::cout << " Hit " << nHit++ << " DetId " << hit->geographicalId().det();
3233  if (hit->geographicalId().det() == DetId::Tracker)
3234  std::cout << " Tracker ";
3235  else if (hit->geographicalId().det() == DetId::Muon)
3236  std::cout << " Muon ";
3237  else
3238  std::cout << " Unknown ";
3239  if (!hit->isValid())
3240  std::cout << " not valid " << std::endl;
3241  else
3242  std::cout << std::endl;
3243  }
3244 }
gather_cfg.cout
cout
Definition: gather_cfg.py:144

◆ debugTrackHit() [2/2]

void GlobalTrackerMuonAlignment::debugTrackHit ( const std::string  title,
reco::TransientTrack &  alongTr 
)

Definition at line 3209 of file GlobalTrackerMuonAlignment.cc.

References gather_cfg::cout, mps_fire::i, DetId::Muon, reco::TransientTrack::recHitsBegin(), reco::TransientTrack::recHitsEnd(), runGCPTkAlMap::title, and DetId::Tracker.

3209  {
3210  std::cout << " ------- " << title << " --------" << std::endl;
3211  int nHit = 1;
3212  for (trackingRecHit_iterator i = alongTr.recHitsBegin(); i != alongTr.recHitsEnd(); i++) {
3213  std::cout << " Hit " << nHit++ << " DetId " << (*i)->geographicalId().det();
3214  if ((*i)->geographicalId().det() == DetId::Tracker)
3215  std::cout << " Tracker ";
3216  else if ((*i)->geographicalId().det() == DetId::Muon)
3217  std::cout << " Muon ";
3218  else
3219  std::cout << " Unknown ";
3220  if (!(*i)->isValid())
3221  std::cout << " not valid " << std::endl;
3222  else
3223  std::cout << std::endl;
3224  }
3225 }
mps_fire.i
i
Definition: mps_fire.py:428
reco::TransientTrack::recHitsEnd
trackingRecHit_iterator recHitsEnd() const
last iterator to RecHits
Definition: TransientTrack.h:126
reco::TransientTrack::recHitsBegin
trackingRecHit_iterator recHitsBegin() const
first iterator to RecHits
Definition: TransientTrack.h:124
edm::OwnVector::const_iterator
Definition: OwnVector.h:41
gather_cfg.cout
cout
Definition: gather_cfg.py:144

◆ debugTrajectory()

void GlobalTrackerMuonAlignment::debugTrajectory ( const std::string  title,
Trajectory &  traj 
)

Definition at line 3299 of file GlobalTrackerMuonAlignment.cc.

References alongMomentum, Trajectory::chiSquared(), gather_cfg::cout, debugTrajectorySOSv(), Trajectory::direction(), Trajectory::firstMeasurement(), Trajectory::foundHits(), Trajectory::isValid(), Trajectory::lastMeasurement(), runGCPTkAlMap::title, and TrajectoryMeasurement::updatedState().

Referenced by muonFitter(), and trackFitter().

3299  {
3300  std::cout << "\n"
3301  << " ...... " << title << " ...... " << std::endl;
3302  if (!traj.isValid()) {
3303  std::cout << " Not valid !!!!!!!! " << std::endl;
3304  return;
3305  }
3306  std::cout << " chi2/Nhit " << traj.chiSquared() << " / " << traj.foundHits();
3307  if (traj.direction() == alongMomentum)
3308  std::cout << " alongMomentum >>>>" << std::endl;
3309  else
3310  std::cout << " oppositeToMomentum <<<<" << std::endl;
3311  this->debugTrajectorySOSv(" firstMeasurementTSOS ", traj.firstMeasurement().updatedState());
3312  //this->debugTrajectorySOSv(" firstMeasurementTSOS ",traj.firstMeasurement().predictedState());
3313  this->debugTrajectorySOSv(" lastMeasurementTSOS ", traj.lastMeasurement().updatedState());
3314  //this->debugTrajectorySOSv(" geom InnermostState", traj.geometricalInnermostState());
3315  std::cout << " . . . . . . . . . . . . . . . . . . . . . . . . . . . . \n" << std::endl;
3316 }
GlobalTrackerMuonAlignment::debugTrajectorySOSv
void debugTrajectorySOSv(const std::string, TrajectoryStateOnSurface)
Definition: GlobalTrackerMuonAlignment.cc:3273
Trajectory::isValid
bool isValid() const
Definition: Trajectory.h:257
Trajectory::chiSquared
float chiSquared() const
Definition: Trajectory.h:241
Trajectory::foundHits
int foundHits() const
Definition: Trajectory.h:206
Trajectory::lastMeasurement
TrajectoryMeasurement const & lastMeasurement() const
Definition: Trajectory.h:150
Trajectory::direction
PropagationDirection const & direction() const
Definition: Trajectory.cc:133
TrajectoryMeasurement::updatedState
TrajectoryStateOnSurface const & updatedState() const
Definition: TrajectoryMeasurement.h:184
Trajectory::firstMeasurement
TrajectoryMeasurement const & firstMeasurement() const
Definition: Trajectory.h:166
gather_cfg.cout
cout
Definition: gather_cfg.py:144

◆ debugTrajectorySOS()

void GlobalTrackerMuonAlignment::debugTrajectorySOS ( const std::string  title,
TrajectoryStateOnSurface &  trajSOS 
)

Definition at line 3247 of file GlobalTrackerMuonAlignment.cc.

References TrajectoryStateOnSurface::cartesianError(), TrajectoryStateOnSurface::charge(), gather_cfg::cout, TrajectoryStateOnSurface::globalMomentum(), TrajectoryStateOnSurface::globalParameters(), TrajectoryStateOnSurface::globalPosition(), TrajectoryStateOnSurface::isValid(), TrajectoryStateOnSurface::localError(), TrajectoryStateOnSurface::localParameters(), PV3DBase< T, PVType, FrameType >::mag(), CartesianTrajectoryError::matrix(), LocalTrajectoryError::matrix(), PV3DBase< T, PVType, FrameType >::perp(), mathSSE::sqrt(), runGCPTkAlMap::title, GlobalTrajectoryParameters::vector(), and LocalTrajectoryParameters::vector().

Referenced by analyzeTrackTrajectory(), muonFitter(), and trackFitter().

3247  {
3248  std::cout << " --- " << title << " --- " << std::endl;
3249  if (!trajSOS.isValid()) {
3250  std::cout << " Not valid !!!! " << std::endl;
3251  return;
3252  }
3253  std::cout << " R |p| " << trajSOS.globalPosition().perp() << " " << trajSOS.globalMomentum().mag() << " charge "
3254  << trajSOS.charge() << std::endl;
3255  std::cout << " x p " << trajSOS.globalParameters().vector()(0) << " " << trajSOS.globalParameters().vector()(1)
3256  << " " << trajSOS.globalParameters().vector()(2) << " " << trajSOS.globalParameters().vector()(3) << " "
3257  << trajSOS.globalParameters().vector()(4) << " " << trajSOS.globalParameters().vector()(5) << std::endl;
3258  std::cout << " +/- " << sqrt(trajSOS.cartesianError().matrix()(0, 0)) << " "
3259  << sqrt(trajSOS.cartesianError().matrix()(1, 1)) << " " << sqrt(trajSOS.cartesianError().matrix()(2, 2))
3260  << " " << sqrt(trajSOS.cartesianError().matrix()(3, 3)) << " "
3261  << sqrt(trajSOS.cartesianError().matrix()(4, 4)) << " " << sqrt(trajSOS.cartesianError().matrix()(5, 5))
3262  << std::endl;
3263  std::cout << " q/p dxy/dz xy " << trajSOS.localParameters().vector()(0) << " "
3264  << trajSOS.localParameters().vector()(1) << " " << trajSOS.localParameters().vector()(2) << " "
3265  << trajSOS.localParameters().vector()(3) << " " << trajSOS.localParameters().vector()(4) << std::endl;
3266  std::cout << " +/- error " << sqrt(trajSOS.localError().matrix()(0, 0)) << " "
3267  << sqrt(trajSOS.localError().matrix()(1, 1)) << " " << sqrt(trajSOS.localError().matrix()(2, 2)) << " "
3268  << sqrt(trajSOS.localError().matrix()(3, 3)) << " " << sqrt(trajSOS.localError().matrix()(4, 4)) << " "
3269  << std::endl;
3270 }
PV3DBase::perp
T perp() const
Definition: PV3DBase.h:69
TrajectoryStateOnSurface::localError
const LocalTrajectoryError & localError() const
Definition: TrajectoryStateOnSurface.h:77
TrajectoryStateOnSurface::globalParameters
const GlobalTrajectoryParameters & globalParameters() const
Definition: TrajectoryStateOnSurface.h:64
TrajectoryStateOnSurface::localParameters
const LocalTrajectoryParameters & localParameters() const
Definition: TrajectoryStateOnSurface.h:73
TrajectoryStateOnSurface::cartesianError
const CartesianTrajectoryError cartesianError() const
Definition: TrajectoryStateOnSurface.h:71
TrajectoryStateOnSurface::globalPosition
GlobalPoint globalPosition() const
Definition: TrajectoryStateOnSurface.h:65
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
LocalTrajectoryParameters::vector
AlgebraicVector5 vector() const
Definition: LocalTrajectoryParameters.h:120
PV3DBase::mag
T mag() const
Definition: PV3DBase.h:64
TrajectoryStateOnSurface::globalMomentum
GlobalVector globalMomentum() const
Definition: TrajectoryStateOnSurface.h:66
LocalTrajectoryError::matrix
const AlgebraicSymMatrix55 & matrix() const
Definition: LocalTrajectoryError.h:60
CartesianTrajectoryError::matrix
const AlgebraicSymMatrix66 & matrix() const
Definition: CartesianTrajectoryError.h:28
gather_cfg.cout
cout
Definition: gather_cfg.py:144

◆ debugTrajectorySOSv()

void GlobalTrackerMuonAlignment::debugTrajectorySOSv ( const std::string  title,
TrajectoryStateOnSurface  trajSOS 
)

Definition at line 3273 of file GlobalTrackerMuonAlignment.cc.

References TrajectoryStateOnSurface::cartesianError(), TrajectoryStateOnSurface::charge(), gather_cfg::cout, TrajectoryStateOnSurface::globalMomentum(), TrajectoryStateOnSurface::globalParameters(), TrajectoryStateOnSurface::globalPosition(), TrajectoryStateOnSurface::isValid(), TrajectoryStateOnSurface::localError(), TrajectoryStateOnSurface::localParameters(), PV3DBase< T, PVType, FrameType >::mag(), CartesianTrajectoryError::matrix(), LocalTrajectoryError::matrix(), PV3DBase< T, PVType, FrameType >::perp(), mathSSE::sqrt(), runGCPTkAlMap::title, GlobalTrajectoryParameters::vector(), and LocalTrajectoryParameters::vector().

Referenced by debugTrajectory(), muonFitter(), and trackFitter().

3273  {
3274  std::cout << " --- " << title << " --- " << std::endl;
3275  if (!trajSOS.isValid()) {
3276  std::cout << " Not valid !!!! " << std::endl;
3277  return;
3278  }
3279  std::cout << " R |p| " << trajSOS.globalPosition().perp() << " " << trajSOS.globalMomentum().mag() << " charge "
3280  << trajSOS.charge() << std::endl;
3281  std::cout << " x p " << trajSOS.globalParameters().vector()(0) << " " << trajSOS.globalParameters().vector()(1)
3282  << " " << trajSOS.globalParameters().vector()(2) << " " << trajSOS.globalParameters().vector()(3) << " "
3283  << trajSOS.globalParameters().vector()(4) << " " << trajSOS.globalParameters().vector()(5) << std::endl;
3284  std::cout << " +/- " << sqrt(trajSOS.cartesianError().matrix()(0, 0)) << " "
3285  << sqrt(trajSOS.cartesianError().matrix()(1, 1)) << " " << sqrt(trajSOS.cartesianError().matrix()(2, 2))
3286  << " " << sqrt(trajSOS.cartesianError().matrix()(3, 3)) << " "
3287  << sqrt(trajSOS.cartesianError().matrix()(4, 4)) << " " << sqrt(trajSOS.cartesianError().matrix()(5, 5))
3288  << std::endl;
3289  std::cout << " q/p dxy/dz xy " << trajSOS.localParameters().vector()(0) << " "
3290  << trajSOS.localParameters().vector()(1) << " " << trajSOS.localParameters().vector()(2) << " "
3291  << trajSOS.localParameters().vector()(3) << " " << trajSOS.localParameters().vector()(4) << std::endl;
3292  std::cout << " +/- error " << sqrt(trajSOS.localError().matrix()(0, 0)) << " "
3293  << sqrt(trajSOS.localError().matrix()(1, 1)) << " " << sqrt(trajSOS.localError().matrix()(2, 2)) << " "
3294  << sqrt(trajSOS.localError().matrix()(3, 3)) << " " << sqrt(trajSOS.localError().matrix()(4, 4)) << " "
3295  << std::endl;
3296 }
PV3DBase::perp
T perp() const
Definition: PV3DBase.h:69
TrajectoryStateOnSurface::localError
const LocalTrajectoryError & localError() const
Definition: TrajectoryStateOnSurface.h:77
TrajectoryStateOnSurface::globalParameters
const GlobalTrajectoryParameters & globalParameters() const
Definition: TrajectoryStateOnSurface.h:64
TrajectoryStateOnSurface::localParameters
const LocalTrajectoryParameters & localParameters() const
Definition: TrajectoryStateOnSurface.h:73
TrajectoryStateOnSurface::cartesianError
const CartesianTrajectoryError cartesianError() const
Definition: TrajectoryStateOnSurface.h:71
TrajectoryStateOnSurface::globalPosition
GlobalPoint globalPosition() const
Definition: TrajectoryStateOnSurface.h:65
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
LocalTrajectoryParameters::vector
AlgebraicVector5 vector() const
Definition: LocalTrajectoryParameters.h:120
PV3DBase::mag
T mag() const
Definition: PV3DBase.h:64
TrajectoryStateOnSurface::globalMomentum
GlobalVector globalMomentum() const
Definition: TrajectoryStateOnSurface.h:66
LocalTrajectoryError::matrix
const AlgebraicSymMatrix55 & matrix() const
Definition: LocalTrajectoryError.h:60
CartesianTrajectoryError::matrix
const AlgebraicSymMatrix66 & matrix() const
Definition: CartesianTrajectoryError.h:28
gather_cfg.cout
cout
Definition: gather_cfg.py:144

◆ endJob()

void GlobalTrackerMuonAlignment::endJob ( )
overrideprivatevirtual

Reimplemented from edm::one::EDAnalyzerBase.

Definition at line 367 of file GlobalTrackerMuonAlignment.cc.

References collectionCosmic, collectionIsolated, gather_cfg::cout, ztail::d, debug_, file, fitHist(), Gfr, Grad, GradL, Hess, HessL, mps_fire::i, Inf, dqmiolumiharvest::j, dqmdumpme::k, MuGlAngle, MuGlShift, MuSelect, N_event, N_track, MillePedeFileConverter_cfg::out, OutGlobalTxt, smuonTags_, mathSSE::sqrt(), txtOutFile_, writeDB_, and writeGlPosRcd().

Referenced by o2olib.O2ORunMgr::executeJob().

367  {
368  bool alarm = false;
369 
370  this->fitHist();
371 
372  AlgebraicVector3 d(0., 0., 0.); // ------------ alignmnet Global Algebraic
373 
374  AlgebraicSymMatrix33 InfI; // inverse it
375  for (int i = 0; i <= 2; i++)
376  for (int j = 0; j <= 2; j++) {
377  if (j < i)
378  continue;
379  InfI(i, j) += Inf(i, j);
380  }
381  bool ierr = !InfI.Invert();
382  if (ierr) {
383  if (alarm)
384  std::cout << " Error inverse Inf matrix !!!!!!!!!!!" << std::endl;
385  }
386 
387  for (int i = 0; i <= 2; i++)
388  for (int k = 0; k <= 2; k++)
389  d(i) -= InfI(i, k) * Gfr(k);
390  // end of Global Algebraic
391 
392  // --------------- alignment Global CLHEP
393  CLHEP::HepVector d3 = CLHEP::solve(Hess, -Grad);
394  int iEr3;
395  CLHEP::HepMatrix Errd3 = Hess.inverse(iEr3);
396  if (iEr3 != 0) {
397  if (alarm)
398  std::cout << " endJob Error inverse Hess matrix !!!!!!!!!!!" << std::endl;
399  }
400  // end of Global CLHEP
401 
402  // ----------------- alignment Local CLHEP
403  CLHEP::HepVector dLI = CLHEP::solve(HessL, -GradL);
404  int iErI;
405  CLHEP::HepMatrix ErrdLI = HessL.inverse(iErI);
406  if (iErI != 0) {
407  if (alarm)
408  std::cout << " endJob Error inverse HessL matrix !!!!!!!!!!!" << std::endl;
409  }
410  // end of Local CLHEP
411 
412  // printout of final parameters
413  std::cout << " ---- " << N_event << " event " << N_track << " tracks " << MuSelect << " ---- " << std::endl;
414  if (collectionIsolated)
415  std::cout << " ALCARECOMuAlCalIsolatedMu " << std::endl;
416  else if (collectionCosmic)
417  std::cout << " ALCARECOMuAlGlobalCosmics " << std::endl;
418  else
419  std::cout << smuonTags_ << std::endl;
420  std::cout << " Similated shifts[cm] " << MuGlShift(1) << " " << MuGlShift(2) << " " << MuGlShift(3) << " "
421  << " angles[rad] " << MuGlAngle(1) << " " << MuGlAngle(2) << " " << MuGlAngle(3) << " " << std::endl;
422  std::cout << " d " << -d << std::endl;
423  ;
424  std::cout << " +- " << sqrt(InfI(0, 0)) << " " << sqrt(InfI(1, 1)) << " " << sqrt(InfI(2, 2)) << std::endl;
425  std::cout << " dG " << d3(1) << " " << d3(2) << " " << d3(3) << " " << d3(4) << " " << d3(5) << " " << d3(6)
426  << std::endl;
427  ;
428  std::cout << " +- " << sqrt(Errd3(1, 1)) << " " << sqrt(Errd3(2, 2)) << " " << sqrt(Errd3(3, 3)) << " "
429  << sqrt(Errd3(4, 4)) << " " << sqrt(Errd3(5, 5)) << " " << sqrt(Errd3(6, 6)) << std::endl;
430  std::cout << " dL " << dLI(1) << " " << dLI(2) << " " << dLI(3) << " " << dLI(4) << " " << dLI(5) << " " << dLI(6)
431  << std::endl;
432  ;
433  std::cout << " +- " << sqrt(ErrdLI(1, 1)) << " " << sqrt(ErrdLI(2, 2)) << " " << sqrt(ErrdLI(3, 3)) << " "
434  << sqrt(ErrdLI(4, 4)) << " " << sqrt(ErrdLI(5, 5)) << " " << sqrt(ErrdLI(6, 6)) << std::endl;
435 
436  // what do we write to DB
437  CLHEP::HepVector vectorToDb(6, 0), vectorErrToDb(6, 0);
438  //vectorToDb = d3;
439  //for(unsigned int i=1; i<=6; i++) vectorErrToDb(i) = sqrt(Errd3(i,i));
440  vectorToDb = -dLI;
441  for (unsigned int i = 1; i <= 6; i++)
442  vectorErrToDb(i) = sqrt(ErrdLI(i, i));
443 
444  // write histograms to root file
445  file->Write();
446  file->Close();
447 
448  // write global parameters to text file
449  OutGlobalTxt.open(txtOutFile_.c_str(), ios::out);
450  if (!OutGlobalTxt.is_open())
451  std::cout << " outglobal.txt is not open !!!!!" << std::endl;
452  else {
453  OutGlobalTxt << vectorToDb(1) << " " << vectorToDb(2) << " " << vectorToDb(3) << " " << vectorToDb(4) << " "
454  << vectorToDb(5) << " " << vectorToDb(6) << " muon Global.\n";
455  OutGlobalTxt << vectorErrToDb(1) << " " << vectorErrToDb(1) << " " << vectorErrToDb(1) << " " << vectorErrToDb(1)
456  << " " << vectorErrToDb(1) << " " << vectorErrToDb(1) << " errors.\n";
457  OutGlobalTxt << N_event << " events are processed.\n";
458 
459  if (collectionIsolated)
460  OutGlobalTxt << "ALCARECOMuAlCalIsolatedMu.\n";
461  else if (collectionCosmic)
462  OutGlobalTxt << " ALCARECOMuAlGlobalCosmics.\n";
463  else
464  OutGlobalTxt << smuonTags_ << ".\n";
465  OutGlobalTxt.close();
466  std::cout << " Write to the file outglobal.txt done " << std::endl;
467  }
468 
469  // write new GlobalPositionRcd to DB
470  if (debug_)
471  std::cout << " writeBD_ " << writeDB_ << std::endl;
472  if (writeDB_)
473  GlobalTrackerMuonAlignment::writeGlPosRcd(vectorToDb);
474 }
mps_fire.i
i
Definition: mps_fire.py:428
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
ztail.d
d
Definition: ztail.py:151
GlobalTrackerMuonAlignment::writeGlPosRcd
void writeGlPosRcd(CLHEP::HepVector &d3)
Definition: GlobalTrackerMuonAlignment.cc:3319
MillePedeFileConverter_cfg.out
out
Definition: MillePedeFileConverter_cfg.py:31
AlgebraicVector3
ROOT::Math::SVector< double, 3 > AlgebraicVector3
Definition: AlgebraicROOTObjects.h:12
AlgebraicSymMatrix33
ROOT::Math::SMatrix< double, 3, 3, ROOT::Math::MatRepSym< double, 3 > > AlgebraicSymMatrix33
Definition: AlgebraicROOTObjects.h:21
gather_cfg.cout
cout
Definition: gather_cfg.py:144

◆ fitHist()

void GlobalTrackerMuonAlignment::fitHist ( )

Definition at line 555 of file GlobalTrackerMuonAlignment.cc.

References histo12, histo13, histo14, histo15, histo16, histo17, histo21, histo22, histo23, histo24, histo25, histo26, histo32, histo33, histo34, histo35, and histo7.

Referenced by endJob().

555  {
556  histo7->Fit("gaus", "Q");
557 
558  histo12->Fit("gaus", "Q");
559  histo13->Fit("gaus", "Q");
560  histo14->Fit("gaus", "Q");
561  histo15->Fit("gaus", "Q");
562  histo16->Fit("gaus", "Q");
563  histo17->Fit("gaus", "Q");
564 
565  histo21->Fit("gaus", "Q");
566  histo22->Fit("gaus", "Q");
567  histo23->Fit("gaus", "Q");
568  histo24->Fit("gaus", "Q");
569  histo25->Fit("gaus", "Q");
570  histo26->Fit("gaus", "Q");
571 
572  histo32->Fit("gaus", "Q");
573  histo33->Fit("gaus", "Q");
574  histo34->Fit("gaus", "Q");
575  histo35->Fit("gaus", "Q");
576 }

◆ gradientGlobal()

void GlobalTrackerMuonAlignment::gradientGlobal ( GlobalVector &  GRt,
GlobalVector &  GPt,
GlobalVector &  GRm,
GlobalVector &  GPm,
GlobalVector &  GNorm,
AlgebraicSymMatrix66 &  GCov 
)

Definition at line 2120 of file GlobalTrackerMuonAlignment.cc.

References A, a, B, CLHEP_dot(), gather_cfg::cout, ztail::d, debug_, dumpMFGeometry_cfg::delta, MillePedeFileConverter_cfg::e, Grad, Hess, mps_fire::i, cuy::ii, info(), dqmiolumiharvest::j, AlCaHLTBitMon_ParallelJobs::p, alignCSCRings::r, alignCSCRings::s, mathSSE::sqrt(), cms::cuda::V, w(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by analyzeTrackTrack(), and analyzeTrackTrajectory().

2125  {
2126  // we search for 6D global correction vector (d, a), where
2127  // 3D vector of shihts d
2128  // 3D vector of rotation angles a
2129 
2130  //bool alarm = true;
2131  bool alarm = false;
2132  bool info = false;
2133 
2134  int Nd = 6; // dimension of vector of alignment pararmeters, d
2135 
2136  //double PtMom = GPt.mag();
2137  CLHEP::HepSymMatrix w(Nd, 0);
2138  for (int i = 1; i <= Nd; i++)
2139  for (int j = 1; j <= Nd; j++) {
2140  if (j <= i)
2141  w(i, j) = GCov(i - 1, j - 1);
2142  //if(i >= 3) w(i,j) /= PtMom;
2143  //if(j >= 3) w(i,j) /= PtMom;
2144  if ((i == j) && (i <= 3) && (GCov(i - 1, j - 1) < 1.e-20))
2145  w(i, j) = 1.e20; // w=0
2146  if (i != j)
2147  w(i, j) = 0.; // use diaginal elements
2148  }
2149 
2150  //GPt /= GPt.mag();
2151  //GPm /= GPm.mag(); // end of transform
2152 
2153  CLHEP::HepVector V(Nd), Rt(3), Pt(3), Rm(3), Pm(3), Norm(3);
2154  Rt(1) = GRt.x();
2155  Rt(2) = GRt.y();
2156  Rt(3) = GRt.z();
2157  Pt(1) = GPt.x();
2158  Pt(2) = GPt.y();
2159  Pt(3) = GPt.z();
2160  Rm(1) = GRm.x();
2161  Rm(2) = GRm.y();
2162  Rm(3) = GRm.z();
2163  Pm(1) = GPm.x();
2164  Pm(2) = GPm.y();
2165  Pm(3) = GPm.z();
2166  Norm(1) = GNorm.x();
2167  Norm(2) = GNorm.y();
2168  Norm(3) = GNorm.z();
2169 
2170  V = dsum(Rm - Rt, Pm - Pt); //std::cout<<" V "<<V.T()<<std::endl;
2171 
2172  //double PmN = CLHEP::dot(Pm, Norm);
2173  double PmN = CLHEP_dot(Pm, Norm);
2174 
2175  CLHEP::HepMatrix Jac(Nd, Nd, 0);
2176  for (int i = 1; i <= 3; i++)
2177  for (int j = 1; j <= 3; j++) {
2178  Jac(i, j) = Pm(i) * Norm(j) / PmN;
2179  if (i == j)
2180  Jac(i, j) -= 1.;
2181  }
2182 
2183  // dp/da
2184  Jac(4, 4) = 0.; // dpx/dax
2185  Jac(5, 4) = -Pm(3); // dpy/dax
2186  Jac(6, 4) = Pm(2); // dpz/dax
2187  Jac(4, 5) = Pm(3); // dpx/day
2188  Jac(5, 5) = 0.; // dpy/day
2189  Jac(6, 5) = -Pm(1); // dpz/day
2190  Jac(4, 6) = -Pm(2); // dpx/daz
2191  Jac(5, 6) = Pm(1); // dpy/daz
2192  Jac(6, 6) = 0.; // dpz/daz
2193 
2194  CLHEP::HepVector dsda(3);
2195  dsda(1) = (Norm(2) * Rm(3) - Norm(3) * Rm(2)) / PmN;
2196  dsda(2) = (Norm(3) * Rm(1) - Norm(1) * Rm(3)) / PmN;
2197  dsda(3) = (Norm(1) * Rm(2) - Norm(2) * Rm(1)) / PmN;
2198 
2199  // dr/da
2200  Jac(1, 4) = Pm(1) * dsda(1); // drx/dax
2201  Jac(2, 4) = -Rm(3) + Pm(2) * dsda(1); // dry/dax
2202  Jac(3, 4) = Rm(2) + Pm(3) * dsda(1); // drz/dax
2203 
2204  Jac(1, 5) = Rm(3) + Pm(1) * dsda(2); // drx/day
2205  Jac(2, 5) = Pm(2) * dsda(2); // dry/day
2206  Jac(3, 5) = -Rm(1) + Pm(3) * dsda(2); // drz/day
2207 
2208  Jac(1, 6) = -Rm(2) + Pm(1) * dsda(3); // drx/daz
2209  Jac(2, 6) = Rm(1) + Pm(2) * dsda(3); // dry/daz
2210  Jac(3, 6) = Pm(3) * dsda(3); // drz/daz
2211 
2212  CLHEP::HepSymMatrix W(Nd, 0);
2213  int ierr;
2214  W = w.inverse(ierr);
2215  if (ierr != 0) {
2216  if (alarm)
2217  std::cout << " gradientGlobal: inversion of matrix w fail " << std::endl;
2218  return;
2219  }
2220 
2221  CLHEP::HepMatrix W_Jac(Nd, Nd, 0);
2222  W_Jac = Jac.T() * W;
2223 
2224  CLHEP::HepVector grad3(Nd);
2225  grad3 = W_Jac * V;
2226 
2227  CLHEP::HepMatrix hess3(Nd, Nd);
2228  hess3 = Jac.T() * W * Jac;
2229  //hess3(4,4) = 1.e-10; hess3(5,5) = 1.e-10; hess3(6,6) = 1.e-10; //????????????????
2230 
2231  Grad += grad3;
2232  Hess += hess3;
2233 
2234  CLHEP::HepVector d3I = CLHEP::solve(Hess, -Grad);
2235  int iEr3I;
2236  CLHEP::HepMatrix Errd3I = Hess.inverse(iEr3I);
2237  if (iEr3I != 0) {
2238  if (alarm)
2239  std::cout << " gradientGlobal error inverse Hess matrix !!!!!!!!!!!" << std::endl;
2240  }
2241 
2242  if (info || debug_) {
2243  std::cout << " dG " << d3I(1) << " " << d3I(2) << " " << d3I(3) << " " << d3I(4) << " " << d3I(5) << " " << d3I(6)
2244  << std::endl;
2245  ;
2246  std::cout << " +- " << sqrt(Errd3I(1, 1)) << " " << sqrt(Errd3I(2, 2)) << " " << sqrt(Errd3I(3, 3)) << " "
2247  << sqrt(Errd3I(4, 4)) << " " << sqrt(Errd3I(5, 5)) << " " << sqrt(Errd3I(6, 6)) << std::endl;
2248  }
2249 
2250 #ifdef CHECK_OF_DERIVATIVES
2251  // -------------------- check of derivatives
2252 
2253  CLHEP::HepVector d(3, 0), a(3, 0);
2254  CLHEP::HepMatrix T(3, 3, 1);
2255 
2256  CLHEP::HepMatrix Ti = T.T();
2257  //double A = CLHEP::dot(Ti*Pm, Norm);
2258  //double B = CLHEP::dot((Rt -Ti*Rm + Ti*d), Norm);
2259  double A = CLHEP_dot(Ti * Pm, Norm);
2260  double B = CLHEP_dot((Rt - Ti * Rm + Ti * d), Norm);
2261  double s0 = B / A;
2262 
2263  CLHEP::HepVector r0(3, 0), p0(3, 0);
2264  r0 = Ti * Rm - Ti * d + s0 * (Ti * Pm) - Rt;
2265  p0 = Ti * Pm - Pt;
2266 
2267  double delta = 0.0001;
2268 
2269  int ii = 3;
2270  d(ii) += delta; // d
2271  //T(2,3) += delta; T(3,2) -= delta; int ii = 1; // a1
2272  //T(3,1) += delta; T(1,3) -= delta; int ii = 2; // a2
2273  //T(1,2) += delta; T(2,1) -= delta; int ii = 3; // a2
2274  Ti = T.T();
2275  //A = CLHEP::dot(Ti*Pm, Norm);
2276  //B = CLHEP::dot((Rt -Ti*Rm + Ti*d), Norm);
2277  A = CLHEP_dot(Ti * Pm, Norm);
2278  B = CLHEP_dot((Rt - Ti * Rm + Ti * d), Norm);
2279  double s = B / A;
2280 
2281  CLHEP::HepVector r(3, 0), p(3, 0);
2282  r = Ti * Rm - Ti * d + s * (Ti * Pm) - Rt;
2283  p = Ti * Pm - Pt;
2284 
2285  std::cout << " s0 s num dsda(" << ii << ") " << s0 << " " << s << " " << (s - s0) / delta << " " << dsda(ii) << endl;
2286  // d(r,p) / d shift
2287  std::cout << " -- An d(r,p)/d(" << ii << ") " << Jac(1, ii) << " " << Jac(2, ii) << " " << Jac(3, ii) << " "
2288  << Jac(4, ii) << " " << Jac(5, ii) << " " << Jac(6, ii) << std::endl;
2289  std::cout << " Nu d(r,p)/d(" << ii << ") " << (r(1) - r0(1)) / delta << " " << (r(2) - r0(2)) / delta << " "
2290  << (r(3) - r0(3)) / delta << " " << (p(1) - p0(1)) / delta << " " << (p(2) - p0(2)) / delta << " "
2291  << (p(3) - p0(3)) / delta << std::endl;
2292  // d(r,p) / d angle
2293  std::cout << " -- An d(r,p)/a(" << ii << ") " << Jac(1, ii + 3) << " " << Jac(2, ii + 3) << " " << Jac(3, ii + 3)
2294  << " " << Jac(4, ii + 3) << " " << Jac(5, ii + 3) << " " << Jac(6, ii + 3) << std::endl;
2295  std::cout << " Nu d(r,p)/a(" << ii << ") " << (r(1) - r0(1)) / delta << " " << (r(2) - r0(2)) / delta << " "
2296  << (r(3) - r0(3)) / delta << " " << (p(1) - p0(1)) / delta << " " << (p(2) - p0(2)) / delta << " "
2297  << (p(3) - p0(3)) / delta << std::endl;
2298  // ----------------------------- end of check
2299 #endif
2300 
2301  return;
2302 } // end gradientGlobal
GlobalTrackerMuonAlignment::CLHEP_dot
double CLHEP_dot(const CLHEP::HepVector &a, const CLHEP::HepVector &b)
Definition: GlobalTrackerMuonAlignment.cc:148
info
static const TGPicture * info(bool iBackgroundIsBlack)
Definition: FWCollectionSummaryWidget.cc:153
mps_fire.i
i
Definition: mps_fire.py:428
B
Definition: APVGainStruct.h:7
w
T w() const
Definition: extBasic3DVector.h:225
PV3DBase::z
T z() const
Definition: PV3DBase.h:61
MillePedeFileConverter_cfg.e
e
Definition: MillePedeFileConverter_cfg.py:37
dumpMFGeometry_cfg.delta
delta
Definition: dumpMFGeometry_cfg.py:25
PV3DBase::x
T x() const
Definition: PV3DBase.h:59
PV3DBase::y
T y() const
Definition: PV3DBase.h:60
cms::cuda::V
uint32_t T const *__restrict__ uint32_t const *__restrict__ int32_t int Histo::index_type cudaStream_t V
Definition: HistoContainer.h:51
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
ztail.d
d
Definition: ztail.py:151
cuy.ii
ii
Definition: cuy.py:589
a
double a
Definition: hdecay.h:119
gather_cfg.cout
cout
Definition: gather_cfg.py:144
A
Definition: APVGainStruct.h:7
T
long double T
Definition: Basic3DVectorLD.h:48
AlCaHLTBitMon_ParallelJobs.p
def p
Definition: AlCaHLTBitMon_ParallelJobs.py:153

◆ gradientGlobalAlg()

void GlobalTrackerMuonAlignment::gradientGlobalAlg ( GlobalVector &  Rt,
GlobalVector &  Pt,
GlobalVector &  Rm,
GlobalVector &  Nl,
AlgebraicSymMatrix66 &  Cm 
)

Definition at line 2029 of file GlobalTrackerMuonAlignment.cc.

References gather_cfg::cout, debug_, HGC3DClusterGenMatchSelector_cfi::dR, MillePedeFileConverter_cfg::e, Gfr, mps_fire::i, Inf, dqmiolumiharvest::j, dqmdumpme::k, mathSSE::sqrt(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by analyzeTrackTrack(), and analyzeTrackTrajectory().

2030  {
2031  // ---------------------------- Calculate Information matrix and Gfree vector
2032  // Information == Hessian , Gfree == gradient of the objective function
2033 
2034  AlgebraicMatrix33 Jac;
2035  AlgebraicVector3 Wi, R_m, R_t, P_t, Norm, dR;
2036 
2037  R_m(0) = Rm.x();
2038  R_m(1) = Rm.y();
2039  R_m(2) = Rm.z();
2040  R_t(0) = Rt.x();
2041  R_t(1) = Rt.y();
2042  R_t(2) = Rt.z();
2043  P_t(0) = Pt.x();
2044  P_t(1) = Pt.y();
2045  P_t(2) = Pt.z();
2046  Norm(0) = Nl.x();
2047  Norm(1) = Nl.y();
2048  Norm(2) = Nl.z();
2049 
2050  for (int i = 0; i <= 2; i++) {
2051  if (Cm(i, i) > 1.e-20)
2052  Wi(i) = 1. / Cm(i, i);
2053  else
2054  Wi(i) = 1.e-10;
2055  dR(i) = R_m(i) - R_t(i);
2056  }
2057 
2058  float PtN = P_t(0) * Norm(0) + P_t(1) * Norm(1) + P_t(2) * Norm(2);
2059 
2060  Jac(0, 0) = 1. - P_t(0) * Norm(0) / PtN;
2061  Jac(0, 1) = -P_t(0) * Norm(1) / PtN;
2062  Jac(0, 2) = -P_t(0) * Norm(2) / PtN;
2063 
2064  Jac(1, 0) = -P_t(1) * Norm(0) / PtN;
2065  Jac(1, 1) = 1. - P_t(1) * Norm(1) / PtN;
2066  Jac(1, 2) = -P_t(1) * Norm(2) / PtN;
2067 
2068  Jac(2, 0) = -P_t(2) * Norm(0) / PtN;
2069  Jac(2, 1) = -P_t(2) * Norm(1) / PtN;
2070  Jac(2, 2) = 1. - P_t(2) * Norm(2) / PtN;
2071 
2072  AlgebraicSymMatrix33 Itr;
2073 
2074  for (int i = 0; i <= 2; i++)
2075  for (int j = 0; j <= 2; j++) {
2076  if (j < i)
2077  continue;
2078  Itr(i, j) = 0.;
2079  //std::cout<<" ij "<<i<<" "<<j<<std::endl;
2080  for (int k = 0; k <= 2; k++) {
2081  Itr(i, j) += Jac(k, i) * Wi(k) * Jac(k, j);
2082  }
2083  }
2084 
2085  for (int i = 0; i <= 2; i++)
2086  for (int j = 0; j <= 2; j++) {
2087  if (j < i)
2088  continue;
2089  Inf(i, j) += Itr(i, j);
2090  }
2091 
2092  AlgebraicVector3 Gtr(0., 0., 0.);
2093  for (int i = 0; i <= 2; i++)
2094  for (int k = 0; k <= 2; k++)
2095  Gtr(i) += dR(k) * Wi(k) * Jac(k, i);
2096  for (int i = 0; i <= 2; i++)
2097  Gfr(i) += Gtr(i);
2098 
2099  if (debug_) {
2100  std::cout << " Wi " << Wi << std::endl;
2101  std::cout << " N " << Norm << std::endl;
2102  std::cout << " P_t " << P_t << std::endl;
2103  std::cout << " (Pt*N) " << PtN << std::endl;
2104  std::cout << " dR " << dR << std::endl;
2105  std::cout << " +/- " << 1. / sqrt(Wi(0)) << " " << 1. / sqrt(Wi(1)) << " " << 1. / sqrt(Wi(2)) << " " << std::endl;
2106  std::cout << " Jacobian dr/ddx " << std::endl;
2107  std::cout << Jac << std::endl;
2108  std::cout << " G-- " << Gtr << std::endl;
2109  std::cout << " Itrack " << std::endl;
2110  std::cout << Itr << std::endl;
2111  std::cout << " Gfr " << Gfr << std::endl;
2112  std::cout << " -- Inf --" << std::endl;
2113  std::cout << Inf << std::endl;
2114  }
2115 
2116  return;
2117 }
mps_fire.i
i
Definition: mps_fire.py:428
AlgebraicMatrix33
ROOT::Math::SMatrix< double, 3, 3, ROOT::Math::MatRepStd< double, 3, 3 > > AlgebraicMatrix33
Definition: AlgebraicROOTObjects.h:41
PV3DBase::z
T z() const
Definition: PV3DBase.h:61
MillePedeFileConverter_cfg.e
e
Definition: MillePedeFileConverter_cfg.py:37
edm::eventsetup::Itr
RecordProviders::iterator Itr
Definition: EventSetupProvider.cc:355
PV3DBase::x
T x() const
Definition: PV3DBase.h:59
PV3DBase::y
T y() const
Definition: PV3DBase.h:60
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
HGC3DClusterGenMatchSelector_cfi.dR
dR
Definition: HGC3DClusterGenMatchSelector_cfi.py:7
AlgebraicVector3
ROOT::Math::SVector< double, 3 > AlgebraicVector3
Definition: AlgebraicROOTObjects.h:12
AlgebraicSymMatrix33
ROOT::Math::SMatrix< double, 3, 3, ROOT::Math::MatRepSym< double, 3 > > AlgebraicSymMatrix33
Definition: AlgebraicROOTObjects.h:21
gather_cfg.cout
cout
Definition: gather_cfg.py:144

◆ gradientLocal()

void GlobalTrackerMuonAlignment::gradientLocal ( GlobalVector &  GRt,
GlobalVector &  GPt,
GlobalVector &  GRm,
GlobalVector &  GPm,
GlobalVector &  GNorm,
CLHEP::HepSymMatrix &  covLoc,
CLHEP::HepMatrix &  rotLoc,
CLHEP::HepVector &  R0,
AlgebraicVector4 &  LPRm 
)

Definition at line 2305 of file GlobalTrackerMuonAlignment.cc.

References CLHEP_dot(), gather_cfg::cout, debug_, dumpMFGeometry_cfg::delta, GradL, HessL, mps_fire::i, cuy::ii, info(), dqmiolumiharvest::j, HLT_2022v11_cff::R0, mathSSE::sqrt(), cms::cuda::V, PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by analyzeTrackTrack(), and analyzeTrackTrajectory().

2313  {
2314  // we search for 6D global correction vector (d, a), where
2315  // 3D vector of shihts d
2316  // 3D vector of rotation angles a
2317 
2318  bool alarm = true;
2319  //bool alarm = false;
2320  bool info = false;
2321 
2322  if (debug_)
2323  std::cout << " gradientLocal " << std::endl;
2324 
2325  /*
2326  const Surface& refSurface = tsosMuon.surface();
2327 
2328  CLHEP::HepMatrix rotLoc (3,3,0);
2329  rotLoc(1,1) = refSurface.rotation().xx();
2330  rotLoc(1,2) = refSurface.rotation().xy();
2331  rotLoc(1,3) = refSurface.rotation().xz();
2332 
2333  rotLoc(2,1) = refSurface.rotation().yx();
2334  rotLoc(2,2) = refSurface.rotation().yy();
2335  rotLoc(2,3) = refSurface.rotation().yz();
2336 
2337  rotLoc(3,1) = refSurface.rotation().zx();
2338  rotLoc(3,2) = refSurface.rotation().zy();
2339  rotLoc(3,3) = refSurface.rotation().zz();
2340  */
2341 
2342  CLHEP::HepVector Rt(3), Pt(3), Rm(3), Pm(3), Norm(3);
2343  Rt(1) = GRt.x();
2344  Rt(2) = GRt.y();
2345  Rt(3) = GRt.z();
2346  Pt(1) = GPt.x();
2347  Pt(2) = GPt.y();
2348  Pt(3) = GPt.z();
2349  Rm(1) = GRm.x();
2350  Rm(2) = GRm.y();
2351  Rm(3) = GRm.z();
2352  Pm(1) = GPm.x();
2353  Pm(2) = GPm.y();
2354  Pm(3) = GPm.z();
2355  Norm(1) = GNorm.x();
2356  Norm(2) = GNorm.y();
2357  Norm(3) = GNorm.z();
2358 
2359  CLHEP::HepVector V(4), Rml(3), Pml(3), Rtl(3), Ptl(3);
2360 
2361  /*
2362  R0(1) = refSurface.position().x();
2363  R0(2) = refSurface.position().y();
2364  R0(3) = refSurface.position().z();
2365  */
2366 
2367  Rml = rotLoc * (Rm - R0);
2368  Rtl = rotLoc * (Rt - R0);
2369  Pml = rotLoc * Pm;
2370  Ptl = rotLoc * Pt;
2371 
2372  V(1) = LPRm(0) - Ptl(1) / Ptl(3);
2373  V(2) = LPRm(1) - Ptl(2) / Ptl(3);
2374  V(3) = LPRm(2) - Rtl(1);
2375  V(4) = LPRm(3) - Rtl(2);
2376 
2377  /*
2378  CLHEP::HepSymMatrix Cov(4,0), W(4,0);
2379  for(int i=1; i<=4; i++)
2380  for(int j=1; j<=i; j++){
2381  Cov(i,j) = (tsosMuon.localError().matrix()
2382  + tsosTrack.localError().matrix())(i,j);
2383  //if(i != j) Cov(i,j) = 0.;
2384  //if((i == j) && ((i==1) || (i==2))) Cov(i,j) = 100.;
2385  //if((i == j) && ((i==3) || (i==4))) Cov(i,j) = 10000.;
2386  }
2387  W = Cov;
2388  */
2389 
2390  CLHEP::HepSymMatrix W = covLoc;
2391 
2392  int ierr;
2393  W.invert(ierr);
2394  if (ierr != 0) {
2395  if (alarm)
2396  std::cout << " gradientLocal: inversion of matrix W fail " << std::endl;
2397  return;
2398  }
2399 
2400  // JacobianCartesianToLocal
2401 
2402  //AlgebraicMatrix56 jacobian // differ from calculation above
2403  //= JacobianCartesianToLocal::JacobianCartesianToLocal
2404  //(refSurface, tsosTrack.localParameters()).jacobian();
2405  //for(int i=1; i<=4; i++) for(int j=1; j<=6; j++){
2406  //int j1 = j - 1; JacToLoc(i,j) = jacobian(i, j1);}
2407 
2408  CLHEP::HepMatrix JacToLoc(4, 6, 0);
2409  for (int i = 1; i <= 2; i++)
2410  for (int j = 1; j <= 3; j++) {
2411  JacToLoc(i, j + 3) = (rotLoc(i, j) - rotLoc(3, j) * Pml(i) / Pml(3)) / Pml(3);
2412  JacToLoc(i + 2, j) = rotLoc(i, j);
2413  }
2414 
2415  // JacobianCorrectionsToCartesian
2416  //double PmN = CLHEP::dot(Pm, Norm);
2417  double PmN = CLHEP_dot(Pm, Norm);
2418 
2419  CLHEP::HepMatrix Jac(6, 6, 0);
2420  for (int i = 1; i <= 3; i++)
2421  for (int j = 1; j <= 3; j++) {
2422  Jac(i, j) = Pm(i) * Norm(j) / PmN;
2423  if (i == j)
2424  Jac(i, j) -= 1.;
2425  }
2426 
2427  // dp/da
2428  Jac(4, 4) = 0.; // dpx/dax
2429  Jac(5, 4) = -Pm(3); // dpy/dax
2430  Jac(6, 4) = Pm(2); // dpz/dax
2431  Jac(4, 5) = Pm(3); // dpx/day
2432  Jac(5, 5) = 0.; // dpy/day
2433  Jac(6, 5) = -Pm(1); // dpz/day
2434  Jac(4, 6) = -Pm(2); // dpx/daz
2435  Jac(5, 6) = Pm(1); // dpy/daz
2436  Jac(6, 6) = 0.; // dpz/daz
2437 
2438  CLHEP::HepVector dsda(3);
2439  dsda(1) = (Norm(2) * Rm(3) - Norm(3) * Rm(2)) / PmN;
2440  dsda(2) = (Norm(3) * Rm(1) - Norm(1) * Rm(3)) / PmN;
2441  dsda(3) = (Norm(1) * Rm(2) - Norm(2) * Rm(1)) / PmN;
2442 
2443  // dr/da
2444  Jac(1, 4) = Pm(1) * dsda(1); // drx/dax
2445  Jac(2, 4) = -Rm(3) + Pm(2) * dsda(1); // dry/dax
2446  Jac(3, 4) = Rm(2) + Pm(3) * dsda(1); // drz/dax
2447 
2448  Jac(1, 5) = Rm(3) + Pm(1) * dsda(2); // drx/day
2449  Jac(2, 5) = Pm(2) * dsda(2); // dry/day
2450  Jac(3, 5) = -Rm(1) + Pm(3) * dsda(2); // drz/day
2451 
2452  Jac(1, 6) = -Rm(2) + Pm(1) * dsda(3); // drx/daz
2453  Jac(2, 6) = Rm(1) + Pm(2) * dsda(3); // dry/daz
2454  Jac(3, 6) = Pm(3) * dsda(3); // drz/daz
2455 
2456  // JacobianCorrectionToLocal
2457  CLHEP::HepMatrix JacCorLoc(4, 6, 0);
2458  JacCorLoc = JacToLoc * Jac;
2459 
2460  // gradient and Hessian
2461  CLHEP::HepMatrix W_Jac(6, 4, 0);
2462  W_Jac = JacCorLoc.T() * W;
2463 
2464  CLHEP::HepVector gradL(6);
2465  gradL = W_Jac * V;
2466 
2467  CLHEP::HepMatrix hessL(6, 6);
2468  hessL = JacCorLoc.T() * W * JacCorLoc;
2469 
2470  GradL += gradL;
2471  HessL += hessL;
2472 
2473  CLHEP::HepVector dLI = CLHEP::solve(HessL, -GradL);
2474  int iErI;
2475  CLHEP::HepMatrix ErrdLI = HessL.inverse(iErI);
2476  if (iErI != 0) {
2477  if (alarm)
2478  std::cout << " gradLocal Error inverse Hess matrix !!!!!!!!!!!" << std::endl;
2479  }
2480 
2481  if (info || debug_) {
2482  std::cout << " dL " << dLI(1) << " " << dLI(2) << " " << dLI(3) << " " << dLI(4) << " " << dLI(5) << " " << dLI(6)
2483  << std::endl;
2484  ;
2485  std::cout << " +- " << sqrt(ErrdLI(1, 1)) << " " << sqrt(ErrdLI(2, 2)) << " " << sqrt(ErrdLI(3, 3)) << " "
2486  << sqrt(ErrdLI(4, 4)) << " " << sqrt(ErrdLI(5, 5)) << " " << sqrt(ErrdLI(6, 6)) << std::endl;
2487  }
2488 
2489  if (debug_) {
2490  //std::cout<<" dV(da3) {"<<-JacCorLoc(1,6)*0.002<<" "<<-JacCorLoc(2,6)*0.002<<" "
2491  // <<-JacCorLoc(3,6)*0.002<<" "<<-JacCorLoc(4,6)*0.002<<"}"<<std::endl;
2492  //std::cout<<" JacCLo {"<<JacCorLoc(1,6)<<" "<<JacCorLoc(2,6)<<" "
2493  // <<JacCorLoc(3,6)<<" "<<JacCorLoc(4,6)<<"}"<<std::endl;
2494  //std::cout<<"Jpx/yx {"<<Jac(4,6)/Pm(3)<<" "<<Jac(5,6)/Pm(3)<<" "
2495  // <<Jac(1,6)<<" "<<Jac(2,6)<<"}"<<std::endl;
2496  //std::cout<<"Jac(,a3){"<<Jac(1,6)<<" "<<Jac(2,6)<<" "<<Jac(3,6)<<" "<<Jac(4,6)<<" "
2497  // <<Jac(5,6)<<" "<<Jac(6,6)<<std::endl;
2498  int i = 5;
2499  if (GNorm.z() > 0.95)
2500  std::cout << " Ecap1 N " << GNorm << std::endl;
2501  else if (GNorm.z() < -0.95)
2502  std::cout << " Ecap2 N " << GNorm << std::endl;
2503  else
2504  std::cout << " Barrel N " << GNorm << std::endl;
2505  std::cout << " JacCLo(i," << i << ") = {" << JacCorLoc(1, i) << " " << JacCorLoc(2, i) << " " << JacCorLoc(3, i)
2506  << " " << JacCorLoc(4, i) << "}" << std::endl;
2507  std::cout << " rotLoc " << rotLoc << std::endl;
2508  std::cout << " position " << R0 << std::endl;
2509  std::cout << " Pm,l " << Pm.T() << " " << Pml(1) / Pml(3) << " " << Pml(2) / Pml(3) << std::endl;
2510  std::cout << " Pt,l " << Pt.T() << " " << Ptl(1) / Ptl(3) << " " << Ptl(2) / Ptl(3) << std::endl;
2511  std::cout << " V " << V.T() << std::endl;
2512  std::cout << " Cov \n" << covLoc << std::endl;
2513  std::cout << " W*Cov " << W * covLoc << std::endl;
2514  //std::cout<<" JacCarToLoc = drldrc \n"<<
2515  //JacobianCartesianToLocal::JacobianCartesianToLocal
2516  //(refSurface, tsosTrack.localParameters()).jacobian()<<std::endl;
2517  std::cout << " JacToLoc " << JacToLoc << std::endl;
2518  }
2519 
2520 #ifdef CHECK_OF_JACOBIAN_CARTESIAN_TO_LOCAL
2521  //---------------------- check of derivatives
2522  CLHEP::HepVector V0(4, 0);
2523  V0(1) = Pml(1) / Pml(3) - Ptl(1) / Ptl(3);
2524  V0(2) = Pml(2) / Pml(3) - Ptl(2) / Ptl(3);
2525  V0(3) = Rml(1) - Rtl(1);
2526  V0(4) = Rml(2) - Rtl(2);
2527  int ii = 3;
2528  float delta = 0.01;
2529  CLHEP::HepVector V1(4, 0);
2530  if (ii <= 3) {
2531  Rm(ii) += delta;
2532  Rml = rotLoc * (Rm - R0);
2533  } else {
2534  Pm(ii - 3) += delta;
2535  Pml = rotLoc * Pm;
2536  }
2537  //if(ii <= 3) {Rt(ii) += delta; Rtl = rotLoc * (Rt - R0);}
2538  //else {Pt(ii-3) += delta; Ptl = rotLoc * Pt;}
2539  V1(1) = Pml(1) / Pml(3) - Ptl(1) / Ptl(3);
2540  V1(2) = Pml(2) / Pml(3) - Ptl(2) / Ptl(3);
2541  V1(3) = Rml(1) - Rtl(1);
2542  V1(4) = Rml(2) - Rtl(2);
2543 
2544  if (GNorm.z() > 0.95)
2545  std::cout << " Ecap1 N " << GNorm << std::endl;
2546  else if (GNorm.z() < -0.95)
2547  std::cout << " Ecap2 N " << GNorm << std::endl;
2548  else
2549  std::cout << " Barrel N " << GNorm << std::endl;
2550  std::cout << " dldc Num (i," << ii << ") " << (V1(1) - V0(1)) / delta << " " << (V1(2) - V0(2)) / delta << " "
2551  << (V1(3) - V0(3)) / delta << " " << (V1(4) - V0(4)) / delta << std::endl;
2552  std::cout << " dldc Ana (i," << ii << ") " << JacToLoc(1, ii) << " " << JacToLoc(2, ii) << " " << JacToLoc(3, ii)
2553  << " " << JacToLoc(4, ii) << std::endl;
2554  //float dtxdpx = (rotLoc(1,1) - rotLoc(3,1)*Pml(1)/Pml(3))/Pml(3);
2555  //float dtydpx = (rotLoc(2,1) - rotLoc(3,2)*Pml(2)/Pml(3))/Pml(3);
2556  //std::cout<<" dtx/dpx dty/ "<<dtxdpx<<" "<<dtydpx<<std::endl;
2557 #endif
2558 
2559  return;
2560 } // end gradientLocal
GlobalTrackerMuonAlignment::CLHEP_dot
double CLHEP_dot(const CLHEP::HepVector &a, const CLHEP::HepVector &b)
Definition: GlobalTrackerMuonAlignment.cc:148
info
static const TGPicture * info(bool iBackgroundIsBlack)
Definition: FWCollectionSummaryWidget.cc:153
mps_fire.i
i
Definition: mps_fire.py:428
PV3DBase::z
T z() const
Definition: PV3DBase.h:61
dumpMFGeometry_cfg.delta
delta
Definition: dumpMFGeometry_cfg.py:25
PV3DBase::x
T x() const
Definition: PV3DBase.h:59
PV3DBase::y
T y() const
Definition: PV3DBase.h:60
cms::cuda::V
uint32_t T const *__restrict__ uint32_t const *__restrict__ int32_t int Histo::index_type cudaStream_t V
Definition: HistoContainer.h:51
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
cuy.ii
ii
Definition: cuy.py:589
HLT_2022v11_cff.R0
R0
Definition: HLT_2022v11_cff.py:6878
gather_cfg.cout
cout
Definition: gather_cfg.py:144

◆ misalignMuon()

void GlobalTrackerMuonAlignment::misalignMuon ( GlobalVector &  GRm,
GlobalVector &  GPm,
GlobalVector &  Nl,
GlobalVector &  Rt,
GlobalVector &  Rm,
GlobalVector &  Pm 
)

Definition at line 2563 of file GlobalTrackerMuonAlignment.cc.

References A, a, B, alignmentValidation::c1, CLHEP_dot(), funct::cos(), gather_cfg::cout, ztail::d, debug_, MuGlAngle, MuGlShift, mkfit::Config::Rl, alignCSCRings::s, funct::sin(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

2564  {
2565  CLHEP::HepVector d(3, 0), a(3, 0), Norm(3), Rm0(3), Pm0(3), Rm1(3), Pm1(3), Rmi(3), Pmi(3);
2566 
2567  d(1) = 0.0;
2568  d(2) = 0.0;
2569  d(3) = 0.0;
2570  a(1) = 0.0000;
2571  a(2) = 0.0000;
2572  a(3) = 0.0000; // zero
2573  //d(1)=0.0; d(2)=0.0; d(3)=0.0; a(1)=0.0020; a(2)=0.0000; a(3)=0.0000; // a1=.002
2574  //d(1)=0.0; d(2)=0.0; d(3)=0.0; a(1)=0.0000; a(2)=0.0020; a(3)=0.0000; // a2=.002
2575  //d(1)=0.0; d(2)=0.0; d(3)=0.0; a(1)=0.0000; a(2)=0.0000; a(3)=0.0020; // a3=.002
2576  //d(1)=1.0; d(2)=2.0; d(3)=3.0; a(1)=0.0000; a(2)=0.0000; a(3)=0.0000; // 12300
2577  //d(1)=0.0; d(2)=0.0; d(3)=0.0; a(1)=0.0000; a(2)=0.0000; a(3)=0.0020; // a3=0.002
2578  //d(1)=10.; d(2)=20.; d(3)=30.; a(1)=0.0100; a(2)=0.0200; a(3)=0.0300; // huge
2579  //d(1)=10.; d(2)=20.; d(3)=30.; a(1)=0.2000; a(2)=0.2500; a(3)=0.3000; // huge angles
2580  //d(1)=0.3; d(2)=0.4; d(3)=0.5; a(1)=0.0005; a(2)=0.0006; a(3)=0.0007; // 0345 0567
2581  //d(1)=0.3; d(2)=0.4; d(3)=0.5; a(1)=0.0001; a(2)=0.0002; a(3)=0.0003; // 0345 0123
2582  //d(1)=0.2; d(2)=0.2; d(3)=0.2; a(1)=0.0002; a(2)=0.0002; a(3)=0.0002; // 0111 0111
2583 
2584  MuGlShift = d;
2585  MuGlAngle = a;
2586 
2587  if ((d(1) == 0.) & (d(2) == 0.) & (d(3) == 0.) & (a(1) == 0.) & (a(2) == 0.) & (a(3) == 0.)) {
2588  Rm = GRm;
2589  Pm = GPm;
2590  if (debug_)
2591  std::cout << " ...... without misalignment " << std::endl;
2592  return;
2593  }
2594 
2595  Rm0(1) = GRm.x();
2596  Rm0(2) = GRm.y();
2597  Rm0(3) = GRm.z();
2598  Pm0(1) = GPm.x();
2599  Pm0(2) = GPm.y();
2600  Pm0(3) = GPm.z();
2601  Norm(1) = Nl.x();
2602  Norm(2) = Nl.y();
2603  Norm(3) = Nl.z();
2604 
2605  CLHEP::HepMatrix T(3, 3, 1);
2606 
2607  //T(1,2) = a(3); T(1,3) = -a(2);
2608  //T(2,1) = -a(3); T(2,3) = a(1);
2609  //T(3,1) = a(2); T(3,2) = -a(1);
2610 
2611  double s1 = std::sin(a(1)), c1 = std::cos(a(1));
2612  double s2 = std::sin(a(2)), c2 = std::cos(a(2));
2613  double s3 = std::sin(a(3)), c3 = std::cos(a(3));
2614  T(1, 1) = c2 * c3;
2615  T(1, 2) = c1 * s3 + s1 * s2 * c3;
2616  T(1, 3) = s1 * s3 - c1 * s2 * c3;
2617  T(2, 1) = -c2 * s3;
2618  T(2, 2) = c1 * c3 - s1 * s2 * s3;
2619  T(2, 3) = s1 * c3 + c1 * s2 * s3;
2620  T(3, 1) = s2;
2621  T(3, 2) = -s1 * c2;
2622  T(3, 3) = c1 * c2;
2623 
2624  //int terr;
2625  //CLHEP::HepMatrix Ti = T.inverse(terr);
2626  CLHEP::HepMatrix Ti = T.T();
2627  CLHEP::HepVector di = -Ti * d;
2628 
2629  CLHEP::HepVector ex0(3, 0), ey0(3, 0), ez0(3, 0), ex(3, 0), ey(3, 0), ez(3, 0);
2630  ex0(1) = 1.;
2631  ey0(2) = 1.;
2632  ez0(3) = 1.;
2633  ex = Ti * ex0;
2634  ey = Ti * ey0;
2635  ez = Ti * ez0;
2636 
2637  CLHEP::HepVector TiN = Ti * Norm;
2638  //double si = CLHEP::dot((Ti*Rm0 - Rm0 + di), TiN) / CLHEP::dot(Pm0, TiN);
2639  //Rm1(1) = CLHEP::dot(ex, Rm0 + si*Pm0 - di);
2640  //Rm1(2) = CLHEP::dot(ey, Rm0 + si*Pm0 - di);
2641  //Rm1(3) = CLHEP::dot(ez, Rm0 + si*Pm0 - di);
2642  double si = CLHEP_dot((Ti * Rm0 - Rm0 + di), TiN) / CLHEP_dot(Pm0, TiN);
2643  Rm1(1) = CLHEP_dot(ex, Rm0 + si * Pm0 - di);
2644  Rm1(2) = CLHEP_dot(ey, Rm0 + si * Pm0 - di);
2645  Rm1(3) = CLHEP_dot(ez, Rm0 + si * Pm0 - di);
2646  Pm1 = T * Pm0;
2647 
2648  Rm = GlobalVector(Rm1(1), Rm1(2), Rm1(3));
2649  //Pm = GlobalVector(CLHEP::dot(Pm0,ex), CLHEP::dot(Pm0, ey), CLHEP::dot(Pm0,ez));// =T*Pm0
2650  Pm = GlobalVector(CLHEP_dot(Pm0, ex), CLHEP_dot(Pm0, ey), CLHEP_dot(Pm0, ez)); // =T*Pm0
2651 
2652  if (debug_) { // ------------- debug tranformation
2653 
2654  std::cout << " ----- Pm " << Pm << std::endl;
2655  std::cout << " T*Pm0 " << Pm1.T() << std::endl;
2656  //std::cout<<" Ti*T*Pm0 "<<(Ti*(T*Pm0)).T()<<std::endl;
2657 
2658  //CLHEP::HepVector Rml = (Rm0 + si*Pm0 - di) + di;
2659  CLHEP::HepVector Rml = Rm1(1) * ex + Rm1(2) * ey + Rm1(3) * ez + di;
2660 
2661  CLHEP::HepVector Rt0(3);
2662  Rt0(1) = Rt.x();
2663  Rt0(2) = Rt.y();
2664  Rt0(3) = Rt.z();
2665 
2666  //double sl = CLHEP::dot(T*(Rt0 - Rml), T*Norm) / CLHEP::dot(Ti * Pm1, Norm);
2667  double sl = CLHEP_dot(T * (Rt0 - Rml), T * Norm) / CLHEP_dot(Ti * Pm1, Norm);
2668  CLHEP::HepVector Rl = Rml + sl * (Ti * Pm1);
2669 
2670  //double A = CLHEP::dot(Ti*Pm1, Norm);
2671  //double B = CLHEP::dot(Rt0, Norm) - CLHEP::dot(Ti*Rm1, Norm)
2672  //+ CLHEP::dot(Ti*d, Norm);
2673  double A = CLHEP_dot(Ti * Pm1, Norm);
2674  double B = CLHEP_dot(Rt0, Norm) - CLHEP_dot(Ti * Rm1, Norm) + CLHEP_dot(Ti * d, Norm);
2675  double s = B / A;
2676  CLHEP::HepVector Dr = Ti * Rm1 - Ti * d + s * (Ti * Pm1) - Rm0;
2677  CLHEP::HepVector Dp = Ti * Pm1 - Pm0;
2678 
2679  CLHEP::HepVector TiR = Ti * Rm0 + di;
2680  CLHEP::HepVector Ri = T * TiR + d;
2681 
2682  std::cout << " --TTi-0 " << (Ri - Rm0).T() << std::endl;
2683  std::cout << " -- Dr " << Dr.T() << endl;
2684  std::cout << " -- Dp " << Dp.T() << endl;
2685  //std::cout<<" ex "<<ex.T()<<endl;
2686  //std::cout<<" ey "<<ey.T()<<endl;
2687  //std::cout<<" ez "<<ez.T()<<endl;
2688  //std::cout<<" ---- T ---- "<<T<<std::endl;
2689  //std::cout<<" ---- Ti ---- "<<Ti<<std::endl;
2690  //std::cout<<" ---- d ---- "<<d.T()<<std::endl;
2691  //std::cout<<" ---- di ---- "<<di.T()<<std::endl;
2692  std::cout << " -- si sl s " << si << " " << sl << " " << s << std::endl;
2693  //std::cout<<" -- si sl "<<si<<" "<<sl<<endl;
2694  //std::cout<<" -- si s "<<si<<" "<<s<<endl;
2695  std::cout << " -- Rml-(Rm0+si*Pm0) " << (Rml - Rm0 - si * Pm0).T() << std::endl;
2696  //std::cout<<" -- Rm0 "<<Rm0.T()<<std::endl;
2697  //std::cout<<" -- Rm1 "<<Rm1.T()<<std::endl;
2698  //std::cout<<" -- Rmi "<<Rmi.T()<<std::endl;
2699  //std::cout<<" --siPm "<<(si*Pm0).T()<<std::endl;
2700  //std::cout<<" --s2Pm "<<(s2*(T * Pm1)).T()<<std::endl;
2701  //std::cout<<" --R1-0 "<<(Rm1-Rm0).T()<<std::endl;
2702  //std::cout<<" --Ri-0 "<<(Rmi-Rm0).T()<<std::endl;
2703  std::cout << " --Rl-0 " << (Rl - Rm0).T() << std::endl;
2704  //std::cout<<" --Pi-0 "<<(Pmi-Pm0).T()<<std::endl;
2705  } // -------- end of debug
2706 
2707  return;
2708 
2709 } // end of misalignMuon
GlobalTrackerMuonAlignment::CLHEP_dot
double CLHEP_dot(const CLHEP::HepVector &a, const CLHEP::HepVector &b)
Definition: GlobalTrackerMuonAlignment.cc:148
B
Definition: APVGainStruct.h:7
PV3DBase::z
T z() const
Definition: PV3DBase.h:61
funct::sin
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
mkfit::Config::Rl
constexpr float Rl[136]
Definition: MaterialEffects.cc:8
PV3DBase::x
T x() const
Definition: PV3DBase.h:59
PV3DBase::y
T y() const
Definition: PV3DBase.h:60
funct::cos
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
ztail.d
d
Definition: ztail.py:151
a
double a
Definition: hdecay.h:119
gather_cfg.cout
cout
Definition: gather_cfg.py:144
A
Definition: APVGainStruct.h:7
T
long double T
Definition: Basic3DVectorLD.h:48
GlobalVector
Global3DVector GlobalVector
Definition: GlobalVector.h:10

◆ misalignMuonL()

void GlobalTrackerMuonAlignment::misalignMuonL ( GlobalVector &  GRm,
GlobalVector &  GPm,
GlobalVector &  Nl,
GlobalVector &  Rt,
GlobalVector &  Rm,
GlobalVector &  Pm,
AlgebraicVector4 &  Vm,
TrajectoryStateOnSurface &  tsosTrack,
TrajectoryStateOnSurface &  tsosMuon 
)

Definition at line 2712 of file GlobalTrackerMuonAlignment.cc.

References A, a, B, alignmentValidation::c1, CLHEP_dot(), funct::cos(), gather_cfg::cout, ztail::d, debug_, cuy::ii, TrajectoryStateOnSurface::localParameters(), MuGlAngle, MuGlShift, GloballyPositioned< T >::position(), LocalTrajectoryParameters::pzSign(), HLT_2022v11_cff::R0, mkfit::Config::Rl, GloballyPositioned< T >::rotation(), alignCSCRings::s, funct::sin(), mathSSE::sqrt(), TrajectoryStateOnSurface::surface(), LocalTrajectoryParameters::vector(), PV3DBase< T, PVType, FrameType >::x(), TkRotation< T >::xx(), TkRotation< T >::xy(), TkRotation< T >::xz(), PV3DBase< T, PVType, FrameType >::y(), TkRotation< T >::yx(), TkRotation< T >::yy(), TkRotation< T >::yz(), PV3DBase< T, PVType, FrameType >::z(), TkRotation< T >::zx(), TkRotation< T >::zy(), and TkRotation< T >::zz().

Referenced by analyzeTrackTrack(), and analyzeTrackTrajectory().

2720  {
2721  CLHEP::HepVector d(3, 0), a(3, 0), Norm(3), Rm0(3), Pm0(3), Rm1(3), Pm1(3), Rmi(3), Pmi(3);
2722 
2723  d(1) = 0.0;
2724  d(2) = 0.0;
2725  d(3) = 0.0;
2726  a(1) = 0.0000;
2727  a(2) = 0.0000;
2728  a(3) = 0.0000; // zero
2729  //d(1)=0.0000001; d(2)=0.0; d(3)=0.0; a(1)=0.0000; a(2)=0.0000; a(3)=0.0000; // nonzero
2730  //d(1)=0.0; d(2)=0.0; d(3)=0.0; a(1)=0.0020; a(2)=0.0000; a(3)=0.0000; // a1=.002
2731  //d(1)=0.0; d(2)=0.0; d(3)=0.0; a(1)=0.0000; a(2)=0.0020; a(3)=0.0000; // a2=.002
2732  //d(1)=0.0; d(2)=0.0; d(3)=0.0; a(1)=0.0000; a(2)=0.0000; a(3)=0.0020; // a3=.002
2733  //d(1)=0.0; d(2)=0.0; d(3)=0.0; a(1)=0.0000; a(2)=0.0000; a(3)=0.0200; // a3=.020
2734  //d(1)=1.0; d(2)=2.0; d(3)=3.0; a(1)=0.0000; a(2)=0.0000; a(3)=0.0000; // 12300
2735  //d(1)=0.0; d(2)=0.0; d(3)=0.0; a(1)=0.0000; a(2)=0.0000; a(3)=0.0020; // a3=0.002
2736  //d(1)=10.; d(2)=20.; d(3)=30.; a(1)=0.0100; a(2)=0.0200; a(3)=0.0300; // huge
2737  //d(1)=10.; d(2)=20.; d(3)=30.; a(1)=0.2000; a(2)=0.2500; a(3)=0.3000; // huge angles
2738  //d(1)=0.3; d(2)=0.4; d(3)=0.5; a(1)=0.0005; a(2)=0.0006; a(3)=0.0007; // 0345 0567
2739  //d(1)=0.3; d(2)=0.4; d(3)=0.5; a(1)=0.0001; a(2)=0.0002; a(3)=0.0003; // 0345 0123
2740  //d(1)=0.1; d(2)=0.2; d(3)=0.3; a(1)=0.0003; a(2)=0.0004; a(3)=0.0005; // 0123 0345
2741 
2742  MuGlShift = d;
2743  MuGlAngle = a;
2744 
2745  if ((d(1) == 0.) & (d(2) == 0.) & (d(3) == 0.) & (a(1) == 0.) & (a(2) == 0.) & (a(3) == 0.)) {
2746  Rm = GRm;
2747  Pm = GPm;
2748  AlgebraicVector4 Vm0;
2749  Vm = AlgebraicVector4(tsosMuon.localParameters().vector()(1),
2750  tsosMuon.localParameters().vector()(2),
2751  tsosMuon.localParameters().vector()(3),
2752  tsosMuon.localParameters().vector()(4));
2753  if (debug_)
2754  std::cout << " ...... without misalignment " << std::endl;
2755  return;
2756  }
2757 
2758  Rm0(1) = GRm.x();
2759  Rm0(2) = GRm.y();
2760  Rm0(3) = GRm.z();
2761  Pm0(1) = GPm.x();
2762  Pm0(2) = GPm.y();
2763  Pm0(3) = GPm.z();
2764  Norm(1) = Nl.x();
2765  Norm(2) = Nl.y();
2766  Norm(3) = Nl.z();
2767 
2768  CLHEP::HepMatrix T(3, 3, 1);
2769 
2770  //T(1,2) = a(3); T(1,3) = -a(2);
2771  //T(2,1) = -a(3); T(2,3) = a(1);
2772  //T(3,1) = a(2); T(3,2) = -a(1);
2773 
2774  double s1 = std::sin(a(1)), c1 = std::cos(a(1));
2775  double s2 = std::sin(a(2)), c2 = std::cos(a(2));
2776  double s3 = std::sin(a(3)), c3 = std::cos(a(3));
2777  T(1, 1) = c2 * c3;
2778  T(1, 2) = c1 * s3 + s1 * s2 * c3;
2779  T(1, 3) = s1 * s3 - c1 * s2 * c3;
2780  T(2, 1) = -c2 * s3;
2781  T(2, 2) = c1 * c3 - s1 * s2 * s3;
2782  T(2, 3) = s1 * c3 + c1 * s2 * s3;
2783  T(3, 1) = s2;
2784  T(3, 2) = -s1 * c2;
2785  T(3, 3) = c1 * c2;
2786 
2787  // Ti, di what we have to apply for misalignment
2788  //int terr;
2789  //CLHEP::HepMatrix Ti = T.inverse(terr);
2790  CLHEP::HepMatrix Ti = T.T();
2791  CLHEP::HepVector di = -Ti * d;
2792 
2793  CLHEP::HepVector ex0(3, 0), ey0(3, 0), ez0(3, 0), ex(3, 0), ey(3, 0), ez(3, 0);
2794  ex0(1) = 1.;
2795  ey0(2) = 1.;
2796  ez0(3) = 1.;
2797  ex = Ti * ex0;
2798  ey = Ti * ey0;
2799  ez = Ti * ez0;
2800 
2801  CLHEP::HepVector TiN = Ti * Norm;
2802  //double si = CLHEP::dot((Ti*Rm0 - Rm0 + di), TiN) / CLHEP::dot(Pm0, TiN);
2803  //Rm1(1) = CLHEP::dot(ex, Rm0 + si*Pm0 - di);
2804  //Rm1(2) = CLHEP::dot(ey, Rm0 + si*Pm0 - di);
2805  //Rm1(3) = CLHEP::dot(ez, Rm0 + si*Pm0 - di);
2806  double si = CLHEP_dot((Ti * Rm0 - Rm0 + di), TiN) / CLHEP_dot(Pm0, TiN);
2807  Rm1(1) = CLHEP_dot(ex, Rm0 + si * Pm0 - di);
2808  Rm1(2) = CLHEP_dot(ey, Rm0 + si * Pm0 - di);
2809  Rm1(3) = CLHEP_dot(ez, Rm0 + si * Pm0 - di);
2810  Pm1 = T * Pm0;
2811 
2812  // Global muon parameters after misalignment
2813  Rm = GlobalVector(Rm1(1), Rm1(2), Rm1(3));
2814  //Pm = GlobalVector(CLHEP::dot(Pm0,ex), CLHEP::dot(Pm0, ey), CLHEP::dot(Pm0,ez));// =T*Pm0
2815  Pm = GlobalVector(CLHEP_dot(Pm0, ex), CLHEP_dot(Pm0, ey), CLHEP_dot(Pm0, ez)); // =T*Pm0
2816 
2817  // Local muon parameters after misalignment
2818  const Surface& refSurface = tsosMuon.surface();
2819  CLHEP::HepMatrix Tl(3, 3, 0);
2820 
2821  Tl(1, 1) = refSurface.rotation().xx();
2822  Tl(1, 2) = refSurface.rotation().xy();
2823  Tl(1, 3) = refSurface.rotation().xz();
2824 
2825  Tl(2, 1) = refSurface.rotation().yx();
2826  Tl(2, 2) = refSurface.rotation().yy();
2827  Tl(2, 3) = refSurface.rotation().yz();
2828 
2829  Tl(3, 1) = refSurface.rotation().zx();
2830  Tl(3, 2) = refSurface.rotation().zy();
2831  Tl(3, 3) = refSurface.rotation().zz();
2832 
2833  CLHEP::HepVector R0(3, 0), newR0(3, 0), newRl(3, 0), newPl(3, 0);
2834  R0(1) = refSurface.position().x();
2835  R0(2) = refSurface.position().y();
2836  R0(3) = refSurface.position().z();
2837 
2838  newRl = Tl * (Rm1 - R0);
2839  newPl = Tl * Pm1;
2840 
2841  Vm(0) = newPl(1) / newPl(3);
2842  Vm(1) = newPl(2) / newPl(3);
2843  Vm(2) = newRl(1);
2844  Vm(3) = newRl(2);
2845 
2846 #ifdef CHECK_DERIVATIVES_LOCAL_VS_ANGLES
2847  float del = 0.0001;
2848  //int ii = 1; T(3,2) +=del; T(2,3) -=del;
2849  int ii = 2;
2850  T(3, 1) -= del;
2851  T(1, 3) += del;
2852  //int ii = 3; T(1,2) -=del; T(2,1) +=del;
2853  AlgebraicVector4 Vm0 = Vm;
2854  CLHEP::HepVector Rm10 = Rm1, Pm10 = Pm1;
2855  ;
2856  CLHEP::HepMatrix newTl = Tl * T;
2857  Ti = T.T();
2858  di = -Ti * d;
2859  ex = Ti * ex0;
2860  ey = Ti * ey0;
2861  ez = Ti * ez0;
2862  TiN = Ti * Norm;
2863  si = CLHEP_dot((Ti * Rm0 - Rm0 + di), TiN) / CLHEP_dot(Pm0, TiN);
2864  Rm1(1) = CLHEP_dot(ex, Rm0 + si * Pm0 - di);
2865  Rm1(2) = CLHEP_dot(ey, Rm0 + si * Pm0 - di);
2866  Rm1(3) = CLHEP_dot(ez, Rm0 + si * Pm0 - di);
2867  Pm1 = T * Pm0;
2868 
2869  newRl = Tl * (Rm1 - R0);
2870  newPl = Tl * Pm1;
2871 
2872  Vm(0) = newPl(1) / newPl(3);
2873  Vm(1) = newPl(2) / newPl(3);
2874  Vm(2) = newRl(1);
2875  Vm(3) = newRl(2);
2876  std::cout << " ========= dVm/da" << ii << " " << (Vm - Vm0) * (1. / del) << std::endl;
2877 #endif
2878 
2879  if (debug_) {
2880  //std::cout<<" dRm/da3 "<<((Rm1-Rm10)*(1./del)).T()<<" "<<((Pm1-Pm10)*(1./del)).T()<<std::endl;
2881  //std::cout<<" Norm "<<Norm.T()<<std::endl;
2882  std::cout << " ex " << (Tl.T() * ex0).T() << std::endl;
2883  std::cout << " ey " << (Tl.T() * ey0).T() << std::endl;
2884  std::cout << " ez " << (Tl.T() * ez0).T() << std::endl;
2885  //std::cpot<<"
2886 
2887  std::cout << " pxyz/p gl 0 " << (Pm0 / Pm0.norm()).T() << std::endl;
2888  std::cout << " pxyz/p loc0 " << (Tl * Pm0 / (Tl * Pm0).norm()).T() << std::endl;
2889  std::cout << " pxyz/p glob " << (Pm1 / Pm1.norm()).T() << std::endl;
2890  std::cout << " pxyz/p loc " << (newPl / newPl.norm()).T() << std::endl;
2891 
2892  //std::cout<<" Rot "<<refSurface.rotation()<<endl;
2893  //std::cout<<" Tl "<<Tl<<endl;
2894  std::cout << " ---- phi g0 g1 l0 l1 " << atan2(Pm0(2), Pm0(1)) << " " << atan2(Pm1(2), Pm1(1)) << " "
2895  << atan2((Tl * Pm0)(2), (Tl * Pm0)(1)) << " " << atan2(newPl(2), newPl(1)) << std::endl;
2896  std::cout << " ---- angl Gl Loc " << atan2(Pm1(2), Pm1(1)) - atan2(Pm0(2), Pm0(1)) << " "
2897  << atan2(newPl(2), newPl(1)) - atan2((Tl * Pm0)(2), (Tl * Pm0)(1)) << " "
2898  << atan2(newPl(3), newPl(2)) - atan2((Tl * Pm0)(3), (Tl * Pm0)(2)) << " "
2899  << atan2(newPl(1), newPl(3)) - atan2((Tl * Pm0)(1), (Tl * Pm0)(3)) << " " << std::endl;
2900  }
2901 
2902  if (debug_) {
2903  CLHEP::HepMatrix newTl(3, 3, 0);
2904  CLHEP::HepVector R0(3, 0), newRl(3, 0), newPl(3, 0);
2905  newTl = Tl * Ti.T();
2906  newR0 = Ti * R0 + di;
2907 
2908  std::cout << " N " << Norm.T() << std::endl;
2909  std::cout << " dtxl yl " << Vm(0) - tsosMuon.localParameters().vector()(1) << " "
2910  << Vm(1) - tsosMuon.localParameters().vector()(2) << std::endl;
2911  std::cout << " dXl dYl " << Vm(2) - tsosMuon.localParameters().vector()(3) << " "
2912  << Vm(3) - tsosMuon.localParameters().vector()(4) << std::endl;
2913  std::cout << " localM " << tsosMuon.localParameters().vector() << std::endl;
2914  std::cout << " Vm " << Vm << std::endl;
2915 
2916  CLHEP::HepVector Rvc(3, 0), Pvc(3, 0), Rvg(3, 0), Pvg(3, 0);
2917  Rvc(1) = Vm(2);
2918  Rvc(2) = Vm(3);
2919  Pvc(3) = tsosMuon.localParameters().pzSign() * Pm0.norm() / sqrt(1 + Vm(0) * Vm(0) + Vm(1) * Vm(1));
2920  Pvc(1) = Pvc(3) * Vm(0);
2921  Pvc(2) = Pvc(3) * Vm(1);
2922 
2923  Rvg = newTl.T() * Rvc + newR0;
2924  Pvg = newTl.T() * Pvc;
2925 
2926  std::cout << " newPl " << newPl.T() << std::endl;
2927  std::cout << " Pvc " << Pvc.T() << std::endl;
2928  std::cout << " Rvg " << Rvg.T() << std::endl;
2929  std::cout << " Rm1 " << Rm1.T() << std::endl;
2930  std::cout << " Pvg " << Pvg.T() << std::endl;
2931  std::cout << " Pm1 " << Pm1.T() << std::endl;
2932  }
2933 
2934  if (debug_) { // ---------- how to transform cartesian from shifted to correct
2935 
2936  std::cout << " ----- Pm " << Pm << std::endl;
2937  std::cout << " T*Pm0 " << Pm1.T() << std::endl;
2938  //std::cout<<" Ti*T*Pm0 "<<(Ti*(T*Pm0)).T()<<std::endl;
2939 
2940  //CLHEP::HepVector Rml = (Rm0 + si*Pm0 - di) + di;
2941  CLHEP::HepVector Rml = Rm1(1) * ex + Rm1(2) * ey + Rm1(3) * ez + di;
2942 
2943  CLHEP::HepVector Rt0(3);
2944  Rt0(1) = Rt.x();
2945  Rt0(2) = Rt.y();
2946  Rt0(3) = Rt.z();
2947 
2948  //double sl = CLHEP::dot(T*(Rt0 - Rml), T*Norm) / CLHEP::dot(Ti * Pm1, Norm);
2949  double sl = CLHEP_dot(T * (Rt0 - Rml), T * Norm) / CLHEP_dot(Ti * Pm1, Norm);
2950  CLHEP::HepVector Rl = Rml + sl * (Ti * Pm1);
2951 
2952  //double A = CLHEP::dot(Ti*Pm1, Norm);
2953  //double B = CLHEP::dot(Rt0, Norm) - CLHEP::dot(Ti*Rm1, Norm)
2954  //+ CLHEP::dot(Ti*d, Norm);
2955  double A = CLHEP_dot(Ti * Pm1, Norm);
2956  double B = CLHEP_dot(Rt0, Norm) - CLHEP_dot(Ti * Rm1, Norm) + CLHEP_dot(Ti * d, Norm);
2957  double s = B / A;
2958  CLHEP::HepVector Dr = Ti * Rm1 - Ti * d + s * (Ti * Pm1) - Rm0;
2959  CLHEP::HepVector Dp = Ti * Pm1 - Pm0;
2960 
2961  CLHEP::HepVector TiR = Ti * Rm0 + di;
2962  CLHEP::HepVector Ri = T * TiR + d;
2963 
2964  std::cout << " --TTi-0 " << (Ri - Rm0).T() << std::endl;
2965  std::cout << " -- Dr " << Dr.T() << endl;
2966  std::cout << " -- Dp " << Dp.T() << endl;
2967  //std::cout<<" ex "<<ex.T()<<endl;
2968  //std::cout<<" ey "<<ey.T()<<endl;
2969  //std::cout<<" ez "<<ez.T()<<endl;
2970  //std::cout<<" ---- T ---- "<<T<<std::endl;
2971  //std::cout<<" ---- Ti ---- "<<Ti<<std::endl;
2972  //std::cout<<" ---- d ---- "<<d.T()<<std::endl;
2973  //std::cout<<" ---- di ---- "<<di.T()<<std::endl;
2974  std::cout << " -- si sl s " << si << " " << sl << " " << s << std::endl;
2975  //std::cout<<" -- si sl "<<si<<" "<<sl<<endl;
2976  //std::cout<<" -- si s "<<si<<" "<<s<<endl;
2977  std::cout << " -- Rml-(Rm0+si*Pm0) " << (Rml - Rm0 - si * Pm0).T() << std::endl;
2978  //std::cout<<" -- Rm0 "<<Rm0.T()<<std::endl;
2979  //std::cout<<" -- Rm1 "<<Rm1.T()<<std::endl;
2980  //std::cout<<" -- Rmi "<<Rmi.T()<<std::endl;
2981  //std::cout<<" --siPm "<<(si*Pm0).T()<<std::endl;
2982  //std::cout<<" --s2Pm "<<(s2*(T * Pm1)).T()<<std::endl;
2983  //std::cout<<" --R1-0 "<<(Rm1-Rm0).T()<<std::endl;
2984  //std::cout<<" --Ri-0 "<<(Rmi-Rm0).T()<<std::endl;
2985  std::cout << " --Rl-0 " << (Rl - Rm0).T() << std::endl;
2986  //std::cout<<" --Pi-0 "<<(Pmi-Pm0).T()<<std::endl;
2987  } // -------- end of debug
2988 
2989  return;
2990 
2991 } // end misalignMuonL
LocalTrajectoryParameters::pzSign
float pzSign() const
Sign of the z-component of the momentum in the local frame.
Definition: LocalTrajectoryParameters.h:146
GlobalTrackerMuonAlignment::CLHEP_dot
double CLHEP_dot(const CLHEP::HepVector &a, const CLHEP::HepVector &b)
Definition: GlobalTrackerMuonAlignment.cc:148
B
Definition: APVGainStruct.h:7
TkRotation::xx
T xx() const
Definition: extTkRotation.h:193
PV3DBase::z
T z() const
Definition: PV3DBase.h:61
TkRotation::xy
T xy() const
Definition: extTkRotation.h:194
funct::sin
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
mkfit::Config::Rl
constexpr float Rl[136]
Definition: MaterialEffects.cc:8
TkRotation::zz
T zz() const
Definition: extTkRotation.h:201
TkRotation::yy
T yy() const
Definition: extTkRotation.h:197
TrajectoryStateOnSurface::localParameters
const LocalTrajectoryParameters & localParameters() const
Definition: TrajectoryStateOnSurface.h:73
TkRotation::yz
T yz() const
Definition: extTkRotation.h:198
TrajectoryStateOnSurface::surface
const SurfaceType & surface() const
Definition: TrajectoryStateOnSurface.h:78
TkRotation::zx
T zx() const
Definition: extTkRotation.h:199
PV3DBase::x
T x() const
Definition: PV3DBase.h:59
PV3DBase::y
T y() const
Definition: PV3DBase.h:60
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
LocalTrajectoryParameters::vector
AlgebraicVector5 vector() const
Definition: LocalTrajectoryParameters.h:120
funct::cos
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
TkRotation::zy
T zy() const
Definition: extTkRotation.h:200
AlgebraicVector4
ROOT::Math::SVector< double, 4 > AlgebraicVector4
Definition: AlgebraicROOTObjects.h:13
ztail.d
d
Definition: ztail.py:151
cuy.ii
ii
Definition: cuy.py:589
GloballyPositioned::position
const PositionType & position() const
Definition: GloballyPositioned.h:36
HLT_2022v11_cff.R0
R0
Definition: HLT_2022v11_cff.py:6878
a
double a
Definition: hdecay.h:119
GloballyPositioned::rotation
const RotationType & rotation() const
Definition: GloballyPositioned.h:38
gather_cfg.cout
cout
Definition: gather_cfg.py:144
A
Definition: APVGainStruct.h:7
T
long double T
Definition: Basic3DVectorLD.h:48
TkRotation::xz
T xz() const
Definition: extTkRotation.h:195
TkRotation::yx
T yx() const
Definition: extTkRotation.h:196
GlobalVector
Global3DVector GlobalVector
Definition: GlobalVector.h:10

◆ muonFitter()

void GlobalTrackerMuonAlignment::muonFitter ( reco::TrackRef  alongTr,
reco::TransientTrack &  alongTTr,
PropagationDirection  direction,
TrajectoryStateOnSurface &  trackFittedTSOS 
)

Definition at line 3109 of file GlobalTrackerMuonAlignment.cc.

References alongMomentum, MuonTransientTrackingRecHitBuilder::build(), edm::OwnVector< T, P >::clear(), gather_cfg::cout, debugTrackHit(), debugTrajectory(), debugTrajectorySOS(), debugTrajectorySOSv(), info(), reco::TransientTrack::innermostMeasurementState(), TrajectoryStateOnSurface::localError(), TrajectoryStateOnSurface::localParameters(), magneticField_, LocalTrajectoryError::matrix(), MuRHBuilder, reco::TransientTrack::outermostMeasurementState(), trajectoryStateTransform::persistentState(), DetId::rawId(), rpcPointValidation_cfi::recHit, reco::TransientTrack::recHitsBegin(), reco::TransientTrack::recHitsEnd(), edm::OwnVector< T, P >::size(), TrajectoryStateOnSurface::surface(), theFitter, theFitterOp, theSmoother, theSmootherOp, and TrajectorySmoother::trajectories().

Referenced by analyzeTrackTrajectory().

3112  {
3113  bool info = false;
3114  bool smooth = false;
3115 
3116  if (info)
3117  std::cout << " ......... start of muonFitter ........" << std::endl;
3118  if (false && info) {
3119  this->debugTrackHit(" Muon track hits ", alongTr);
3120  this->debugTrackHit(" Muon TransientTrack hits ", alongTTr);
3121  }
3122 
3123  edm::OwnVector<TrackingRecHit> recHit;
3124  DetId trackDetId = DetId(alongTr->innerDetId());
3125  for (auto const& hit : alongTr->recHits())
3126  recHit.push_back(hit->clone());
3127  if (direction != alongMomentum) {
3128  edm::OwnVector<TrackingRecHit> recHitAlong = recHit;
3129  recHit.clear();
3130  for (unsigned int ihit = recHitAlong.size() - 1; ihit + 1 > 0; ihit--) {
3131  recHit.push_back(recHitAlong[ihit]);
3132  }
3133  recHitAlong.clear();
3134  }
3135  if (info)
3136  std::cout << " ... Own recHit.size() DetId==Muon " << recHit.size() << " " << trackDetId.rawId() << std::endl;
3137 
3138  MuonTransientTrackingRecHitBuilder::ConstRecHitContainer recHitMu =
3139  MuRHBuilder->build(alongTTr.recHitsBegin(), alongTTr.recHitsEnd());
3140  if (info)
3141  std::cout << " ...along.... recHitMu.size() " << recHitMu.size() << std::endl;
3142  if (direction != alongMomentum) {
3143  MuonTransientTrackingRecHitBuilder::ConstRecHitContainer recHitMuAlong = recHitMu;
3144  recHitMu.clear();
3145  for (unsigned int ihit = recHitMuAlong.size() - 1; ihit + 1 > 0; ihit--) {
3146  recHitMu.push_back(recHitMuAlong[ihit]);
3147  }
3148  recHitMuAlong.clear();
3149  }
3150  if (info)
3151  std::cout << " ...opposite... recHitMu.size() " << recHitMu.size() << std::endl;
3152 
3153  TrajectoryStateOnSurface firstTSOS;
3154  if (direction == alongMomentum)
3155  firstTSOS = alongTTr.innermostMeasurementState();
3156  else
3157  firstTSOS = alongTTr.outermostMeasurementState();
3158 
3159  AlgebraicSymMatrix55 CovLoc;
3160  CovLoc(0, 0) = 1. * firstTSOS.localError().matrix()(0, 0);
3161  CovLoc(1, 1) = 10. * firstTSOS.localError().matrix()(1, 1);
3162  CovLoc(2, 2) = 10. * firstTSOS.localError().matrix()(2, 2);
3163  CovLoc(3, 3) = 100. * firstTSOS.localError().matrix()(3, 3);
3164  CovLoc(4, 4) = 100. * firstTSOS.localError().matrix()(4, 4);
3165  TrajectoryStateOnSurface initialTSOS(
3166  firstTSOS.localParameters(), LocalTrajectoryError(CovLoc), firstTSOS.surface(), &*magneticField_);
3167 
3168  PTrajectoryStateOnDet PTraj = trajectoryStateTransform::persistentState(initialTSOS, trackDetId.rawId());
3169  //const TrajectorySeed seedT(*PTraj, recHit, alongMomentum);
3170  const TrajectorySeed seedT(PTraj, recHit, direction);
3171 
3172  std::vector<Trajectory> trajVec;
3173  if (direction == alongMomentum)
3174  trajVec = theFitter->fit(seedT, recHitMu, initialTSOS);
3175  else
3176  trajVec = theFitterOp->fit(seedT, recHitMu, initialTSOS);
3177 
3178  if (info) {
3179  this->debugTrajectorySOSv(" innermostTSOS of TransientTrack", alongTTr.innermostMeasurementState());
3180  this->debugTrajectorySOSv(" outermostTSOS of TransientTrack", alongTTr.outermostMeasurementState());
3181  this->debugTrajectorySOS(" initialTSOS for theFitter ", initialTSOS);
3182  std::cout << " . . . . . trajVec.size() " << trajVec.size() << std::endl;
3183  if (!trajVec.empty())
3184  this->debugTrajectory(" theFitter trajectory ", trajVec[0]);
3185  }
3186 
3187  if (!smooth) {
3188  if (!trajVec.empty())
3189  trackFittedTSOS = trajVec[0].lastMeasurement().updatedState();
3190  } else {
3191  std::vector<Trajectory> trajSVec;
3192  if (!trajVec.empty()) {
3193  if (direction == alongMomentum)
3194  trajSVec = theSmoother->trajectories(*(trajVec.begin()));
3195  else
3196  trajSVec = theSmootherOp->trajectories(*(trajVec.begin()));
3197  }
3198  if (info)
3199  std::cout << " . . . . trajSVec.size() " << trajSVec.size() << std::endl;
3200  if (!trajSVec.empty())
3201  this->debugTrajectorySOSv("smoothed geom InnermostState", trajSVec[0].geometricalInnermostState());
3202  if (!trajSVec.empty())
3203  trackFittedTSOS = trajSVec[0].geometricalInnermostState();
3204  }
3205 
3206 } // end of muonFitter
rpcPointValidation_cfi.recHit
recHit
Definition: rpcPointValidation_cfi.py:7
GlobalTrackerMuonAlignment::debugTrajectorySOSv
void debugTrajectorySOSv(const std::string, TrajectoryStateOnSurface)
Definition: GlobalTrackerMuonAlignment.cc:3273
info
static const TGPicture * info(bool iBackgroundIsBlack)
Definition: FWCollectionSummaryWidget.cc:153
TrajectorySmoother::trajectories
virtual TrajectoryContainer trajectories(const Trajectory &traj) const
Definition: TrajectorySmoother.h:18
TrajectoryStateOnSurface::localError
const LocalTrajectoryError & localError() const
Definition: TrajectoryStateOnSurface.h:77
reco::TransientTrack::recHitsEnd
trackingRecHit_iterator recHitsEnd() const
last iterator to RecHits
Definition: TransientTrack.h:126
MuonTransientTrackingRecHitBuilder::ConstRecHitContainer
TransientTrackingRecHit::ConstRecHitContainer ConstRecHitContainer
Definition: MuonTransientTrackingRecHitBuilder.h:15
trajectoryStateTransform::persistentState
PTrajectoryStateOnDet persistentState(const TrajectoryStateOnSurface &ts, unsigned int detid)
Definition: TrajectoryStateTransform.cc:14
TrajectoryStateOnSurface::localParameters
const LocalTrajectoryParameters & localParameters() const
Definition: TrajectoryStateOnSurface.h:73
reco::TransientTrack::recHitsBegin
trackingRecHit_iterator recHitsBegin() const
first iterator to RecHits
Definition: TransientTrack.h:124
TrajectoryStateOnSurface::surface
const SurfaceType & surface() const
Definition: TrajectoryStateOnSurface.h:78
MuonTransientTrackingRecHitBuilder::build
RecHitPointer build(const TrackingRecHit *p, edm::ESHandle< GlobalTrackingGeometry > trackingGeometry) const
Call the MuonTransientTrackingRecHit::specificBuild.
Definition: MuonTransientTrackingRecHitBuilder.cc:20
GlobalTrackerMuonAlignment::debugTrackHit
void debugTrackHit(const std::string, reco::TrackRef)
Definition: GlobalTrackerMuonAlignment.cc:3228
reco::TransientTrack::outermostMeasurementState
TrajectoryStateOnSurface outermostMeasurementState() const
Definition: TransientTrack.h:84
edm::OwnVector::size
size_type size() const
Definition: OwnVector.h:300
edm::OwnVector::clear
void clear()
Definition: OwnVector.h:481
reco::TransientTrack::innermostMeasurementState
TrajectoryStateOnSurface innermostMeasurementState() const
Definition: TransientTrack.h:86
DetId
Definition: DetId.h:17
DetId::rawId
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:57
AlgebraicSymMatrix55
ROOT::Math::SMatrix< double, 5, 5, ROOT::Math::MatRepSym< double, 5 > > AlgebraicSymMatrix55
Definition: AlgebraicROOTObjects.h:23
GlobalTrackerMuonAlignment::debugTrajectory
void debugTrajectory(const std::string, Trajectory &)
Definition: GlobalTrackerMuonAlignment.cc:3299
LocalTrajectoryError::matrix
const AlgebraicSymMatrix55 & matrix() const
Definition: LocalTrajectoryError.h:60
gather_cfg.cout
cout
Definition: gather_cfg.py:144
GlobalTrackerMuonAlignment::MuRHBuilder
MuonTransientTrackingRecHitBuilder * MuRHBuilder
Definition: GlobalTrackerMuonAlignment.cc:197
GlobalTrackerMuonAlignment::debugTrajectorySOS
void debugTrajectorySOS(const std::string, TrajectoryStateOnSurface &)
Definition: GlobalTrackerMuonAlignment.cc:3247

◆ trackFitter()

void GlobalTrackerMuonAlignment::trackFitter ( reco::TrackRef  alongTr,
reco::TransientTrack &  alongTTr,
PropagationDirection  direction,
TrajectoryStateOnSurface &  trackFittedTSOS 
)

Definition at line 2994 of file GlobalTrackerMuonAlignment.cc.

References alongMomentum, TransientTrackingRecHitBuilder::build(), edm::OwnVector< T, P >::clear(), gather_cfg::cout, debugTrackHit(), debugTrajectory(), debugTrajectorySOS(), debugTrajectorySOSv(), info(), reco::TransientTrack::innermostMeasurementState(), TrajectoryStateOnSurface::localError(), TrajectoryStateOnSurface::localParameters(), magneticField_, LocalTrajectoryError::matrix(), reco::TransientTrack::outermostMeasurementState(), trajectoryStateTransform::persistentState(), DetId::rawId(), rpcPointValidation_cfi::recHit, edm::OwnVector< T, P >::size(), TrajectoryStateOnSurface::surface(), theFitter, theFitterOp, theSmoother, theSmootherOp, TrajectorySmoother::trajectories(), and TTRHBuilder.

Referenced by analyzeTrackTrajectory().

2997  {
2998  bool info = false;
2999  bool smooth = false;
3000 
3001  if (info)
3002  std::cout << " ......... start of trackFitter ......... " << std::endl;
3003  if (false && info) {
3004  this->debugTrackHit(" Tracker track hits ", alongTr);
3005  this->debugTrackHit(" Tracker TransientTrack hits ", alongTTr);
3006  }
3007 
3008  edm::OwnVector<TrackingRecHit> recHit;
3009  DetId trackDetId = DetId(alongTr->innerDetId());
3010  for (auto const& hit : alongTr->recHits())
3011  recHit.push_back(hit->clone());
3012  if (direction != alongMomentum) {
3013  edm::OwnVector<TrackingRecHit> recHitAlong = recHit;
3014  recHit.clear();
3015  for (unsigned int ihit = recHitAlong.size() - 1; ihit + 1 > 0; ihit--) {
3016  recHit.push_back(recHitAlong[ihit]);
3017  }
3018  recHitAlong.clear();
3019  }
3020  if (info)
3021  std::cout << " ... Own recHit.size() " << recHit.size() << " " << trackDetId.rawId() << std::endl;
3022 
3023  //MuonTransientTrackingRecHitBuilder::ConstRecHitContainer recHitTrack;
3024  TransientTrackingRecHit::RecHitContainer recHitTrack;
3025  for (auto const& hit : alongTr->recHits())
3026  recHitTrack.push_back(TTRHBuilder->build(hit));
3027 
3028  if (info)
3029  std::cout << " ... recHitTrack.size() " << recHit.size() << " " << trackDetId.rawId() << " ======> recHitMu "
3030  << std::endl;
3031 
3032  MuonTransientTrackingRecHitBuilder::ConstRecHitContainer recHitMu = recHitTrack;
3033  /*
3034  MuonTransientTrackingRecHitBuilder::ConstRecHitContainer recHitMu =
3035  MuRHBuilder->build(alongTTr.recHits().begin(), alongTTr.recHits().end());
3036  */
3037  if (info)
3038  std::cout << " ...along.... recHitMu.size() " << recHitMu.size() << std::endl;
3039  if (direction != alongMomentum) {
3040  MuonTransientTrackingRecHitBuilder::ConstRecHitContainer recHitMuAlong = recHitMu;
3041  recHitMu.clear();
3042  for (unsigned int ihit = recHitMuAlong.size() - 1; ihit + 1 > 0; ihit--) {
3043  recHitMu.push_back(recHitMuAlong[ihit]);
3044  }
3045  recHitMuAlong.clear();
3046  }
3047  if (info)
3048  std::cout << " ...opposite... recHitMu.size() " << recHitMu.size() << std::endl;
3049 
3050  TrajectoryStateOnSurface firstTSOS;
3051  if (direction == alongMomentum)
3052  firstTSOS = alongTTr.innermostMeasurementState();
3053  else
3054  firstTSOS = alongTTr.outermostMeasurementState();
3055 
3056  AlgebraicSymMatrix55 CovLoc;
3057  CovLoc(0, 0) = 1. * firstTSOS.localError().matrix()(0, 0);
3058  CovLoc(1, 1) = 10. * firstTSOS.localError().matrix()(1, 1);
3059  CovLoc(2, 2) = 10. * firstTSOS.localError().matrix()(2, 2);
3060  CovLoc(3, 3) = 100. * firstTSOS.localError().matrix()(3, 3);
3061  CovLoc(4, 4) = 100. * firstTSOS.localError().matrix()(4, 4);
3062  TrajectoryStateOnSurface initialTSOS(
3063  firstTSOS.localParameters(), LocalTrajectoryError(CovLoc), firstTSOS.surface(), &*magneticField_);
3064 
3065  PTrajectoryStateOnDet PTraj = trajectoryStateTransform::persistentState(initialTSOS, trackDetId.rawId());
3066  //const TrajectorySeed seedT(*PTraj, recHit, alongMomentum);
3067  const TrajectorySeed seedT(PTraj, recHit, direction);
3068 
3069  std::vector<Trajectory> trajVec;
3070  if (direction == alongMomentum)
3071  trajVec = theFitter->fit(seedT, recHitMu, initialTSOS);
3072  else
3073  trajVec = theFitterOp->fit(seedT, recHitMu, initialTSOS);
3074 
3075  if (info) {
3076  this->debugTrajectorySOSv(" innermostTSOS of TransientTrack", alongTTr.innermostMeasurementState());
3077  this->debugTrajectorySOSv(" outermostTSOS of TransientTrack", alongTTr.outermostMeasurementState());
3078  this->debugTrajectorySOS(" initialTSOS for theFitter ", initialTSOS);
3079  std::cout << " . . . . . trajVec.size() " << trajVec.size() << std::endl;
3080  if (!trajVec.empty())
3081  this->debugTrajectory(" theFitter trajectory ", trajVec[0]);
3082  }
3083 
3084  if (!smooth) {
3085  if (!trajVec.empty())
3086  trackFittedTSOS = trajVec[0].lastMeasurement().updatedState();
3087  } else {
3088  std::vector<Trajectory> trajSVec;
3089  if (!trajVec.empty()) {
3090  if (direction == alongMomentum)
3091  trajSVec = theSmoother->trajectories(*(trajVec.begin()));
3092  else
3093  trajSVec = theSmootherOp->trajectories(*(trajVec.begin()));
3094  }
3095  if (info)
3096  std::cout << " . . . . trajSVec.size() " << trajSVec.size() << std::endl;
3097  if (!trajSVec.empty())
3098  this->debugTrajectorySOSv("smoothed geom InnermostState", trajSVec[0].geometricalInnermostState());
3099  if (!trajSVec.empty())
3100  trackFittedTSOS = trajSVec[0].geometricalInnermostState();
3101  }
3102 
3103  if (info)
3104  this->debugTrajectorySOSv(" trackFittedTSOS ", trackFittedTSOS);
3105 
3106 } // end of trackFitter
rpcPointValidation_cfi.recHit
recHit
Definition: rpcPointValidation_cfi.py:7
GlobalTrackerMuonAlignment::debugTrajectorySOSv
void debugTrajectorySOSv(const std::string, TrajectoryStateOnSurface)
Definition: GlobalTrackerMuonAlignment.cc:3273
info
static const TGPicture * info(bool iBackgroundIsBlack)
Definition: FWCollectionSummaryWidget.cc:153
TrajectorySmoother::trajectories
virtual TrajectoryContainer trajectories(const Trajectory &traj) const
Definition: TrajectorySmoother.h:18
TrajectoryStateOnSurface::localError
const LocalTrajectoryError & localError() const
Definition: TrajectoryStateOnSurface.h:77
MuonTransientTrackingRecHitBuilder::ConstRecHitContainer
TransientTrackingRecHit::ConstRecHitContainer ConstRecHitContainer
Definition: MuonTransientTrackingRecHitBuilder.h:15
trajectoryStateTransform::persistentState
PTrajectoryStateOnDet persistentState(const TrajectoryStateOnSurface &ts, unsigned int detid)
Definition: TrajectoryStateTransform.cc:14
TrackingRecHit::RecHitContainer
std::vector< ConstRecHitPointer > RecHitContainer
Definition: TrackingRecHit.h:31
TrajectoryStateOnSurface::localParameters
const LocalTrajectoryParameters & localParameters() const
Definition: TrajectoryStateOnSurface.h:73
TrajectoryStateOnSurface::surface
const SurfaceType & surface() const
Definition: TrajectoryStateOnSurface.h:78
GlobalTrackerMuonAlignment::debugTrackHit
void debugTrackHit(const std::string, reco::TrackRef)
Definition: GlobalTrackerMuonAlignment.cc:3228
reco::TransientTrack::outermostMeasurementState
TrajectoryStateOnSurface outermostMeasurementState() const
Definition: TransientTrack.h:84
TransientTrackingRecHitBuilder::build
virtual RecHitPointer build(const TrackingRecHit *p) const =0
build a tracking rechit from an existing rechit
edm::OwnVector::size
size_type size() const
Definition: OwnVector.h:300
edm::OwnVector::clear
void clear()
Definition: OwnVector.h:481
reco::TransientTrack::innermostMeasurementState
TrajectoryStateOnSurface innermostMeasurementState() const
Definition: TransientTrack.h:86
DetId
Definition: DetId.h:17
DetId::rawId
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:57
GlobalTrackerMuonAlignment::TTRHBuilder
const TransientTrackingRecHitBuilder * TTRHBuilder
Definition: GlobalTrackerMuonAlignment.cc:198
AlgebraicSymMatrix55
ROOT::Math::SMatrix< double, 5, 5, ROOT::Math::MatRepSym< double, 5 > > AlgebraicSymMatrix55
Definition: AlgebraicROOTObjects.h:23
GlobalTrackerMuonAlignment::debugTrajectory
void debugTrajectory(const std::string, Trajectory &)
Definition: GlobalTrackerMuonAlignment.cc:3299
LocalTrajectoryError::matrix
const AlgebraicSymMatrix55 & matrix() const
Definition: LocalTrajectoryError.h:60
gather_cfg.cout
cout
Definition: gather_cfg.py:144
GlobalTrackerMuonAlignment::debugTrajectorySOS
void debugTrajectorySOS(const std::string, TrajectoryStateOnSurface &)
Definition: GlobalTrackerMuonAlignment.cc:3247

◆ writeGlPosRcd()

void GlobalTrackerMuonAlignment::writeGlPosRcd ( CLHEP::HepVector &  d3)

Definition at line 3319 of file GlobalTrackerMuonAlignment.cc.

References alignmentValidation::c1, DetId::Calo, L1TowerCalibrationProducer_cfi::calo, funct::cos(), gather_cfg::cout, cond::service::PoolDBOutputService::currentTime(), debug_, bsc_activity_cfg::ecal, DetId::Ecal, DetId::Hcal, hltEgammaHLTExtra_cfi::hcal, edm::Service< T >::isAvailable(), iteratorEcalRcd, iteratorHcalRcd, iteratorMuonRcd, iteratorTrackerRcd, DetId::Muon, HLT_2022v11_cff::muon, DetId::rawId(), funct::sin(), DetId::Tracker, PbPb_ZMuSkimMuonDPG_cff::tracker, and cond::service::PoolDBOutputService::writeOneIOV().

Referenced by endJob().

3319  {
3320  //paramVec(1) = 0.1; paramVec(2) = 0.2; paramVec(3) = 0.3; //0123
3321  //paramVec(4) = 0.0001; paramVec(5) = 0.0002; paramVec(6) = 0.0003;
3322  //paramVec(1) = 0.3; paramVec(2) = 0.4; paramVec(3) = 0.5; //03450123
3323  //paramVec(4) = 0.0001; paramVec(5) = 0.0002; paramVec(6) = 0.0003;
3324  //paramVec(1) = 2.; paramVec(2) = 2.; paramVec(3) = 2.; //222
3325  //paramVec(4) = 0.02; paramVec(5) = 0.02; paramVec(6) = 0.02;
3326  //paramVec(1) = -10.; paramVec(2) = -20.; paramVec(3) = -30.; //102030
3327  //paramVec(4) = -0.1; paramVec(5) = -0.2; paramVec(6) = -0.3;
3328  //paramVec(1) = 0.; paramVec(2) = 0.; paramVec(3) = 1.; //1z
3329  //paramVec(4) = 0.0000; paramVec(5) = 0.0000; paramVec(6) = 0.0000;
3330 
3331  std::cout << " paramVector " << paramVec.T() << std::endl;
3332 
3333  CLHEP::Hep3Vector colX, colY, colZ;
3334 
3335 #ifdef NOT_EXACT_ROTATION_MATRIX
3336  colX = CLHEP::Hep3Vector(1., -paramVec(6), paramVec(5));
3337  colY = CLHEP::Hep3Vector(paramVec(6), 1., -paramVec(4));
3338  colZ = CLHEP::Hep3Vector(-paramVec(5), paramVec(4), 1.);
3339 #else
3340  double s1 = std::sin(paramVec(4)), c1 = std::cos(paramVec(4));
3341  double s2 = std::sin(paramVec(5)), c2 = std::cos(paramVec(5));
3342  double s3 = std::sin(paramVec(6)), c3 = std::cos(paramVec(6));
3343  colX = CLHEP::Hep3Vector(c2 * c3, -c2 * s3, s2);
3344  colY = CLHEP::Hep3Vector(c1 * s3 + s1 * s2 * c3, c1 * c3 - s1 * s2 * s3, -s1 * c2);
3345  colZ = CLHEP::Hep3Vector(s1 * s3 - c1 * s2 * c3, s1 * c3 + c1 * s2 * s3, c1 * c2);
3346 #endif
3347 
3348  CLHEP::HepVector param0(6, 0);
3349 
3350  Alignments globalPositions{};
3351 
3352  //Tracker
3353  AlignTransform tracker(
3354  iteratorTrackerRcd->translation(), iteratorTrackerRcd->rotation(), DetId(DetId::Tracker).rawId());
3355  // Muon
3356  CLHEP::Hep3Vector posMuGlRcd = iteratorMuonRcd->translation();
3357  CLHEP::HepRotation rotMuGlRcd = iteratorMuonRcd->rotation();
3358  CLHEP::HepEulerAngles angMuGlRcd = iteratorMuonRcd->rotation().eulerAngles();
3359  if (debug_)
3360  std::cout << " Old muon Rcd Euler angles " << angMuGlRcd.phi() << " " << angMuGlRcd.theta() << " "
3361  << angMuGlRcd.psi() << std::endl;
3362  AlignTransform muon;
3363  if ((angMuGlRcd.phi() == 0.) && (angMuGlRcd.theta() == 0.) && (angMuGlRcd.psi() == 0.) && (posMuGlRcd.x() == 0.) &&
3364  (posMuGlRcd.y() == 0.) && (posMuGlRcd.z() == 0.)) {
3365  std::cout << " New muon parameters are stored in Rcd" << std::endl;
3366 
3367  AlignTransform muonNew(AlignTransform::Translation(paramVec(1), paramVec(2), paramVec(3)),
3368  AlignTransform::Rotation(colX, colY, colZ),
3369  DetId(DetId::Muon).rawId());
3370  muon = muonNew;
3371  } else if ((paramVec(1) == 0.) && (paramVec(2) == 0.) && (paramVec(3) == 0.) && (paramVec(4) == 0.) &&
3372  (paramVec(5) == 0.) && (paramVec(6) == 0.)) {
3373  std::cout << " Old muon parameters are stored in Rcd" << std::endl;
3374 
3375  AlignTransform muonNew(iteratorMuonRcd->translation(), iteratorMuonRcd->rotation(), DetId(DetId::Muon).rawId());
3376  muon = muonNew;
3377  } else {
3378  std::cout << " New + Old muon parameters are stored in Rcd" << std::endl;
3379  CLHEP::Hep3Vector posMuGlRcdThis = CLHEP::Hep3Vector(paramVec(1), paramVec(2), paramVec(3));
3380  CLHEP::HepRotation rotMuGlRcdThis = CLHEP::HepRotation(colX, colY, colZ);
3381  CLHEP::Hep3Vector posMuGlRcdNew = rotMuGlRcdThis * posMuGlRcd + posMuGlRcdThis;
3382  CLHEP::HepRotation rotMuGlRcdNew = rotMuGlRcdThis * rotMuGlRcd;
3383 
3384  AlignTransform muonNew(posMuGlRcdNew, rotMuGlRcdNew, DetId(DetId::Muon).rawId());
3385  muon = muonNew;
3386  }
3387 
3388  // Ecal
3389  AlignTransform ecal(iteratorEcalRcd->translation(), iteratorEcalRcd->rotation(), DetId(DetId::Ecal).rawId());
3390  // Hcal
3391  AlignTransform hcal(iteratorHcalRcd->translation(), iteratorHcalRcd->rotation(), DetId(DetId::Hcal).rawId());
3392  // Calo
3393  AlignTransform calo(AlignTransform::Translation(param0(1), param0(2), param0(3)),
3394  AlignTransform::EulerAngles(param0(4), param0(5), param0(6)),
3395  DetId(DetId::Calo).rawId());
3396 
3397  std::cout << "Tracker (" << tracker.rawId() << ") at " << tracker.translation() << " "
3398  << tracker.rotation().eulerAngles() << std::endl;
3399  std::cout << tracker.rotation() << std::endl;
3400 
3401  std::cout << "Muon (" << muon.rawId() << ") at " << muon.translation() << " " << muon.rotation().eulerAngles()
3402  << std::endl;
3403  std::cout << " rotations angles around x,y,z "
3404  << " ( " << -muon.rotation().zy() << " " << muon.rotation().zx() << " " << -muon.rotation().yx() << " )"
3405  << std::endl;
3406  std::cout << muon.rotation() << std::endl;
3407 
3408  std::cout << "Ecal (" << ecal.rawId() << ") at " << ecal.translation() << " " << ecal.rotation().eulerAngles()
3409  << std::endl;
3410  std::cout << ecal.rotation() << std::endl;
3411 
3412  std::cout << "Hcal (" << hcal.rawId() << ") at " << hcal.translation() << " " << hcal.rotation().eulerAngles()
3413  << std::endl;
3414  std::cout << hcal.rotation() << std::endl;
3415 
3416  std::cout << "Calo (" << calo.rawId() << ") at " << calo.translation() << " " << calo.rotation().eulerAngles()
3417  << std::endl;
3418  std::cout << calo.rotation() << std::endl;
3419 
3420  globalPositions.m_align.push_back(tracker);
3421  globalPositions.m_align.push_back(muon);
3422  globalPositions.m_align.push_back(ecal);
3423  globalPositions.m_align.push_back(hcal);
3424  globalPositions.m_align.push_back(calo);
3425 
3426  std::cout << "Uploading to the database..." << std::endl;
3427 
3428  edm::Service<cond::service::PoolDBOutputService> poolDbService;
3429 
3430  if (!poolDbService.isAvailable())
3431  throw cms::Exception("NotAvailable") << "PoolDBOutputService not available";
3432 
3433  // if (poolDbService->isNewTagRequest("GlobalPositionRcd")) {
3434  // poolDbService->createOneIOV<Alignments>(globalPositions, poolDbService->endOfTime(), "GlobalPositionRcd");
3435  // } else {
3436  // poolDbService->appendOneIOV<Alignments>(globalPositions, poolDbService->currentTime(), "GlobalPositionRcd");
3437  // }
3438  poolDbService->writeOneIOV<Alignments>(globalPositions, poolDbService->currentTime(), "GlobalPositionRcd");
3439  std::cout << "done!" << std::endl;
3440 
3441  return;
3442 }
L1TowerCalibrationProducer_cfi.calo
calo
Definition: L1TowerCalibrationProducer_cfi.py:59
HLT_2022v11_cff.muon
muon
Definition: HLT_2022v11_cff.py:9159
funct::sin
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
AlignTransform::Translation
CLHEP::Hep3Vector Translation
Definition: AlignTransform.h:18
GlobalTrackerMuonAlignment::iteratorHcalRcd
std::vector< AlignTransform >::const_iterator iteratorHcalRcd
Definition: GlobalTrackerMuonAlignment.cc:217
AlignTransform::EulerAngles
CLHEP::HepEulerAngles EulerAngles
Definition: AlignTransform.h:17
PbPb_ZMuSkimMuonDPG_cff.tracker
tracker
Definition: PbPb_ZMuSkimMuonDPG_cff.py:60
bsc_activity_cfg.ecal
ecal
Definition: bsc_activity_cfg.py:25
funct::cos
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
GlobalTrackerMuonAlignment::iteratorMuonRcd
std::vector< AlignTransform >::const_iterator iteratorMuonRcd
Definition: GlobalTrackerMuonAlignment.cc:215
cond::service::PoolDBOutputService::writeOneIOV
Hash writeOneIOV(const T &payload, Time_t time, const std::string &recordName)
Definition: PoolDBOutputService.h:62
GlobalTrackerMuonAlignment::iteratorEcalRcd
std::vector< AlignTransform >::const_iterator iteratorEcalRcd
Definition: GlobalTrackerMuonAlignment.cc:216
DetId
Definition: DetId.h:17
DetId::rawId
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:57
GlobalTrackerMuonAlignment::iteratorTrackerRcd
std::vector< AlignTransform >::const_iterator iteratorTrackerRcd
Definition: GlobalTrackerMuonAlignment.cc:214
edm::Service::isAvailable
bool isAvailable() const
Definition: Service.h:40
gather_cfg.cout
cout
Definition: gather_cfg.py:144
hltEgammaHLTExtra_cfi.hcal
hcal
Definition: hltEgammaHLTExtra_cfi.py:42
calo
Definition: Common.h:9
AlignTransform::Rotation
CLHEP::HepRotation Rotation
Definition: AlignTransform.h:19

Member Data Documentation

◆ collectionCosmic

bool GlobalTrackerMuonAlignment::collectionCosmic
private

Definition at line 171 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), beginJob(), and endJob().

◆ collectionIsolated

bool GlobalTrackerMuonAlignment::collectionIsolated
private

Definition at line 172 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), beginJob(), and endJob().

◆ cosmicMuonMode_

bool GlobalTrackerMuonAlignment::cosmicMuonMode_
private

Definition at line 169 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), beginJob(), and GlobalTrackerMuonAlignment().

◆ debug_

bool GlobalTrackerMuonAlignment::debug_
private

Definition at line 178 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), endJob(), gradientGlobal(), gradientGlobalAlg(), gradientLocal(), misalignMuon(), misalignMuonL(), and writeGlPosRcd().

◆ defineFitter

bool GlobalTrackerMuonAlignment::defineFitter
private

Definition at line 196 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrajectory().

◆ file

TFile* GlobalTrackerMuonAlignment::file
private

Definition at line 226 of file GlobalTrackerMuonAlignment.cc.

Referenced by beginJob(), python.rootplot.utilities.RootFile::cd(), endJob(), ztee.GZipLog::finish(), ztee.GZipLog::flush_block(), python.rootplot.utilities.RootFile::get(), tfile.TFileService::stop(), core.VertexHistograms.VertexHistograms::write(), core.TreeAnalyzerNumpy.TreeAnalyzerNumpy::write(), and ztee.GZipLog::write_block().

◆ Gfr

AlgebraicVector3 GlobalTrackerMuonAlignment::Gfr
private

Definition at line 201 of file GlobalTrackerMuonAlignment.cc.

Referenced by beginJob(), endJob(), and gradientGlobalAlg().

◆ globalPositionRcd_

const Alignments* GlobalTrackerMuonAlignment::globalPositionRcd_
private

Definition at line 190 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), and analyzeTrackTrajectory().

◆ gmuonTags_

edm::InputTag GlobalTrackerMuonAlignment::gmuonTags_
private

Definition at line 166 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), and analyzeTrackTrajectory().

◆ Grad

CLHEP::HepVector GlobalTrackerMuonAlignment::Grad
private

Definition at line 204 of file GlobalTrackerMuonAlignment.cc.

Referenced by beginJob(), endJob(), and gradientGlobal().

◆ GradL

CLHEP::HepVector GlobalTrackerMuonAlignment::GradL
private

Definition at line 207 of file GlobalTrackerMuonAlignment.cc.

Referenced by beginJob(), endJob(), and gradientLocal().

◆ Hess

CLHEP::HepMatrix GlobalTrackerMuonAlignment::Hess
private

Definition at line 205 of file GlobalTrackerMuonAlignment.cc.

Referenced by beginJob(), endJob(), and gradientGlobal().

◆ HessL

CLHEP::HepMatrix GlobalTrackerMuonAlignment::HessL
private

Definition at line 208 of file GlobalTrackerMuonAlignment.cc.

Referenced by beginJob(), endJob(), and gradientLocal().

◆ histo

TH1F* GlobalTrackerMuonAlignment::histo
private

Definition at line 227 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and bookHist().

◆ histo101

TH2F* GlobalTrackerMuonAlignment::histo101
private

Definition at line 262 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and bookHist().

◆ histo102

TH2F* GlobalTrackerMuonAlignment::histo102
private

Definition at line 263 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and bookHist().

◆ histo11

TH1F* GlobalTrackerMuonAlignment::histo11
private

Definition at line 236 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and bookHist().

◆ histo12

TH1F* GlobalTrackerMuonAlignment::histo12
private

Definition at line 237 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), bookHist(), and fitHist().

◆ histo13

TH1F* GlobalTrackerMuonAlignment::histo13
private

Definition at line 238 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), bookHist(), and fitHist().

◆ histo14

TH1F* GlobalTrackerMuonAlignment::histo14
private

Definition at line 239 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), bookHist(), and fitHist().

◆ histo15

TH1F* GlobalTrackerMuonAlignment::histo15
private

Definition at line 240 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), bookHist(), and fitHist().

◆ histo16

TH1F* GlobalTrackerMuonAlignment::histo16
private

Definition at line 241 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), bookHist(), and fitHist().

◆ histo17

TH1F* GlobalTrackerMuonAlignment::histo17
private

Definition at line 242 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), bookHist(), and fitHist().

◆ histo18

TH1F* GlobalTrackerMuonAlignment::histo18
private

Definition at line 243 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and bookHist().

◆ histo19

TH1F* GlobalTrackerMuonAlignment::histo19
private

Definition at line 244 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and bookHist().

◆ histo2

TH1F* GlobalTrackerMuonAlignment::histo2
private

Definition at line 228 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and bookHist().

◆ histo20

TH1F* GlobalTrackerMuonAlignment::histo20
private

Definition at line 245 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and bookHist().

◆ histo21

TH1F* GlobalTrackerMuonAlignment::histo21
private

Definition at line 246 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), bookHist(), and fitHist().

◆ histo22

TH1F* GlobalTrackerMuonAlignment::histo22
private

Definition at line 247 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), bookHist(), and fitHist().

◆ histo23

TH1F* GlobalTrackerMuonAlignment::histo23
private

Definition at line 248 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), bookHist(), and fitHist().

◆ histo24

TH1F* GlobalTrackerMuonAlignment::histo24
private

Definition at line 249 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), bookHist(), and fitHist().

◆ histo25

TH1F* GlobalTrackerMuonAlignment::histo25
private

Definition at line 250 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), bookHist(), and fitHist().

◆ histo26

TH1F* GlobalTrackerMuonAlignment::histo26
private

Definition at line 251 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), bookHist(), and fitHist().

◆ histo27

TH1F* GlobalTrackerMuonAlignment::histo27
private

Definition at line 252 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and bookHist().

◆ histo28

TH1F* GlobalTrackerMuonAlignment::histo28
private

Definition at line 253 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and bookHist().

◆ histo29

TH1F* GlobalTrackerMuonAlignment::histo29
private

Definition at line 254 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and bookHist().

◆ histo3

TH1F* GlobalTrackerMuonAlignment::histo3
private

Definition at line 229 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and bookHist().

◆ histo30

TH1F* GlobalTrackerMuonAlignment::histo30
private

Definition at line 255 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and bookHist().

◆ histo31

TH1F* GlobalTrackerMuonAlignment::histo31
private

Definition at line 256 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and bookHist().

◆ histo32

TH1F* GlobalTrackerMuonAlignment::histo32
private

Definition at line 257 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), bookHist(), and fitHist().

◆ histo33

TH1F* GlobalTrackerMuonAlignment::histo33
private

Definition at line 258 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), bookHist(), and fitHist().

◆ histo34

TH1F* GlobalTrackerMuonAlignment::histo34
private

Definition at line 259 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), bookHist(), and fitHist().

◆ histo35

TH1F* GlobalTrackerMuonAlignment::histo35
private

Definition at line 260 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), bookHist(), and fitHist().

◆ histo4

TH1F* GlobalTrackerMuonAlignment::histo4
private

Definition at line 230 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and bookHist().

◆ histo5

TH1F* GlobalTrackerMuonAlignment::histo5
private

Definition at line 231 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and bookHist().

◆ histo6

TH1F* GlobalTrackerMuonAlignment::histo6
private

Definition at line 232 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and bookHist().

◆ histo7

TH1F* GlobalTrackerMuonAlignment::histo7
private

Definition at line 233 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), bookHist(), and fitHist().

◆ histo8

TH1F* GlobalTrackerMuonAlignment::histo8
private

Definition at line 234 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and bookHist().

◆ Inf

AlgebraicSymMatrix33 GlobalTrackerMuonAlignment::Inf
private

Definition at line 202 of file GlobalTrackerMuonAlignment.cc.

Referenced by beginJob(), endJob(), and gradientGlobalAlg().

◆ isolatedMuonMode_

bool GlobalTrackerMuonAlignment::isolatedMuonMode_
private

Definition at line 170 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), beginJob(), and GlobalTrackerMuonAlignment().

◆ iteratorEcalRcd

std::vector<AlignTransform>::const_iterator GlobalTrackerMuonAlignment::iteratorEcalRcd
private

Definition at line 216 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and writeGlPosRcd().

◆ iteratorHcalRcd

std::vector<AlignTransform>::const_iterator GlobalTrackerMuonAlignment::iteratorHcalRcd
private

Definition at line 217 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and writeGlPosRcd().

◆ iteratorMuonRcd

std::vector<AlignTransform>::const_iterator GlobalTrackerMuonAlignment::iteratorMuonRcd
private

Definition at line 215 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and writeGlPosRcd().

◆ iteratorTrackerRcd

std::vector<AlignTransform>::const_iterator GlobalTrackerMuonAlignment::iteratorTrackerRcd
private

Definition at line 214 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and writeGlPosRcd().

◆ m_globalPosToken

edm::ESGetToken<Alignments, GlobalPositionRcd> GlobalTrackerMuonAlignment::m_globalPosToken
private

Definition at line 160 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), and analyzeTrackTrajectory().

◆ m_MagFieldToken

edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> GlobalTrackerMuonAlignment::m_MagFieldToken
private

Definition at line 159 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), and analyzeTrackTrajectory().

◆ m_propToken

edm::ESGetToken<Propagator, TrackingComponentsRecord> GlobalTrackerMuonAlignment::m_propToken
private

Definition at line 161 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), and analyzeTrackTrajectory().

◆ m_TkGeometryToken

edm::ESGetToken<GlobalTrackingGeometry, GlobalTrackingGeometryRecord> GlobalTrackerMuonAlignment::m_TkGeometryToken
private

Definition at line 158 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), and analyzeTrackTrajectory().

◆ m_ttrhBuilderToken

edm::ESGetToken<TransientTrackingRecHitBuilder, TransientRecHitRecord> GlobalTrackerMuonAlignment::m_ttrhBuilderToken
private

Definition at line 162 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrajectory().

◆ magneticField_

const MagneticField* GlobalTrackerMuonAlignment::magneticField_
private

Definition at line 185 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), muonFitter(), and trackFitter().

◆ MuGlAngle

CLHEP::HepVector GlobalTrackerMuonAlignment::MuGlAngle
private

Definition at line 220 of file GlobalTrackerMuonAlignment.cc.

Referenced by endJob(), misalignMuon(), and misalignMuonL().

◆ MuGlShift

CLHEP::HepVector GlobalTrackerMuonAlignment::MuGlShift
private

Definition at line 219 of file GlobalTrackerMuonAlignment.cc.

Referenced by endJob(), misalignMuon(), and misalignMuonL().

◆ muonTags_

edm::InputTag GlobalTrackerMuonAlignment::muonTags_
private

Definition at line 165 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), and analyzeTrackTrajectory().

◆ MuRHBuilder

MuonTransientTrackingRecHitBuilder* GlobalTrackerMuonAlignment::MuRHBuilder
private

Definition at line 197 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrajectory(), and muonFitter().

◆ MuSelect

std::string GlobalTrackerMuonAlignment::MuSelect
private

Definition at line 222 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and endJob().

◆ N_event

int GlobalTrackerMuonAlignment::N_event
private

Definition at line 210 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), beginJob(), and endJob().

◆ N_track

int GlobalTrackerMuonAlignment::N_track
private

Definition at line 211 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), beginJob(), and endJob().

◆ OutGlobalTxt

ofstream GlobalTrackerMuonAlignment::OutGlobalTxt
private

Definition at line 224 of file GlobalTrackerMuonAlignment.cc.

Referenced by endJob().

◆ propagator_

string GlobalTrackerMuonAlignment::propagator_
private

Definition at line 168 of file GlobalTrackerMuonAlignment.cc.

◆ refitMuon_

bool GlobalTrackerMuonAlignment::refitMuon_
private

Definition at line 173 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrajectory().

◆ refitTrack_

bool GlobalTrackerMuonAlignment::refitTrack_
private

Definition at line 174 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrajectory().

◆ rootOutFile_

string GlobalTrackerMuonAlignment::rootOutFile_
private

Definition at line 175 of file GlobalTrackerMuonAlignment.cc.

Referenced by beginJob().

◆ smuonTags_

edm::InputTag GlobalTrackerMuonAlignment::smuonTags_
private

Definition at line 167 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), analyzeTrackTrajectory(), and endJob().

◆ theFitter

KFTrajectoryFitter* GlobalTrackerMuonAlignment::theFitter
private

Definition at line 192 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrajectory(), muonFitter(), and trackFitter().

◆ theFitterOp

KFTrajectoryFitter* GlobalTrackerMuonAlignment::theFitterOp
private

Definition at line 194 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrajectory(), muonFitter(), and trackFitter().

◆ theSmoother

KFTrajectorySmoother* GlobalTrackerMuonAlignment::theSmoother
private

Definition at line 193 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrajectory(), muonFitter(), and trackFitter().

◆ theSmootherOp

KFTrajectorySmoother* GlobalTrackerMuonAlignment::theSmootherOp
private

Definition at line 195 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrajectory(), muonFitter(), and trackFitter().

◆ theTrackingGeometry

edm::ESHandle<GlobalTrackingGeometry> GlobalTrackerMuonAlignment::theTrackingGeometry
private

Definition at line 181 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrajectory().

◆ trackingGeometry_

const GlobalTrackingGeometry* GlobalTrackerMuonAlignment::trackingGeometry_
private

Definition at line 182 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), and analyzeTrackTrajectory().

◆ trackTags_

edm::InputTag GlobalTrackerMuonAlignment::trackTags_
private

Definition at line 164 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), and analyzeTrackTrajectory().

◆ TTRHBuilder

const TransientTrackingRecHitBuilder* GlobalTrackerMuonAlignment::TTRHBuilder
private

Definition at line 198 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrajectory(), and trackFitter().

◆ txtOutFile_

string GlobalTrackerMuonAlignment::txtOutFile_
private

Definition at line 176 of file GlobalTrackerMuonAlignment.cc.

Referenced by endJob().

◆ watchGlobalPositionRcd_

edm::ESWatcher<GlobalPositionRcd> GlobalTrackerMuonAlignment::watchGlobalPositionRcd_
private

Definition at line 189 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), and analyzeTrackTrajectory().

◆ watchMagneticFieldRecord_

edm::ESWatcher<IdealMagneticFieldRecord> GlobalTrackerMuonAlignment::watchMagneticFieldRecord_
private

Definition at line 184 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), and analyzeTrackTrajectory().

◆ watchTrackingComponentsRecord_

edm::ESWatcher<TrackingComponentsRecord> GlobalTrackerMuonAlignment::watchTrackingComponentsRecord_
private

Definition at line 187 of file GlobalTrackerMuonAlignment.cc.

◆ watchTrackingGeometry_

edm::ESWatcher<GlobalTrackingGeometryRecord> GlobalTrackerMuonAlignment::watchTrackingGeometry_
private

Definition at line 180 of file GlobalTrackerMuonAlignment.cc.

Referenced by analyzeTrackTrack(), and analyzeTrackTrajectory().

◆ writeDB_

bool GlobalTrackerMuonAlignment::writeDB_
private

Definition at line 177 of file GlobalTrackerMuonAlignment.cc.

Referenced by endJob().

Generated for CMSSW Reference Manual by doxygen 1.8.14