CMS 3D CMS Logo

List of all members | Public Member Functions | Private Member Functions | Static Private Member Functions | Private Attributes | Static Private Attributes
TkLasBeamFitter Class Reference

#include <Alignment/LaserAlignment/plugins/TkLasBeamFitter.cc>

Inheritance diagram for TkLasBeamFitter:
edm::one::EDProducer< edm::EndRunProducer > edm::one::EDProducerBase edm::ProducerBase edm::EDConsumerBase edm::ProductRegistryHelper

Public Member Functions

void endRunProduce (edm::Run &run, const edm::EventSetup &setup) override
 
void produce (edm::Event &event, const edm::EventSetup &setup) override
 
 TkLasBeamFitter (const edm::ParameterSet &config)
 
 ~TkLasBeamFitter () override
 
- Public Member Functions inherited from edm::one::EDProducer< edm::EndRunProducer >
 EDProducer ()=default
 
 EDProducer (const EDProducer &)=delete
 
SerialTaskQueueglobalLuminosityBlocksQueue () final
 
SerialTaskQueueglobalRunsQueue () final
 
bool hasAbilityToProduceInBeginLumis () const final
 
bool hasAbilityToProduceInBeginProcessBlocks () const final
 
bool hasAbilityToProduceInBeginRuns () const final
 
bool hasAbilityToProduceInEndLumis () const final
 
bool hasAbilityToProduceInEndProcessBlocks () const final
 
bool hasAbilityToProduceInEndRuns () const final
 
const EDProduceroperator= (const EDProducer &)=delete
 
bool wantsGlobalLuminosityBlocks () const final
 
bool wantsGlobalRuns () const final
 
bool wantsInputProcessBlocks () const final
 
bool wantsProcessBlocks () const final
 
- Public Member Functions inherited from edm::one::EDProducerBase
 EDProducerBase ()
 
ModuleDescription const & moduleDescription () const
 
bool wantsStreamLuminosityBlocks () const
 
bool wantsStreamRuns () const
 
 ~EDProducerBase () override
 
- Public Member Functions inherited from edm::ProducerBase
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
std::vector< edm::ProductResolverIndex > const & indiciesForPutProducts (BranchType iBranchType) const
 
 ProducerBase ()
 
std::vector< edm::ProductResolverIndex > const & putTokenIndexToProductResolverIndex () const
 
std::vector< bool > const & recordProvenanceList () const
 
void registerProducts (ProducerBase *, ProductRegistry *, ModuleDescription const &)
 
std::function< void(BranchDescription const &)> registrationCallback () const
 used by the fwk to register list of products More...
 
void resolvePutIndicies (BranchType iBranchType, ModuleToResolverIndicies const &iIndicies, std::string const &moduleLabel)
 
TypeLabelList const & typeLabelList () const
 used by the fwk to register the list of products of this module More...
 
 ~ProducerBase () noexcept(false) override
 
- Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfoconsumesInfo () 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
 
ESResolverIndex const * esGetTokenIndices (edm::Transition iTrans) const
 
std::vector< ESResolverIndex > 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
 
EDConsumerBaseoperator= (EDConsumerBase &&)=default
 
bool registeredToConsume (ProductResolverIndex, bool, 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::ESRecordsToProductResolverIndices const &)
 
virtual ~EDConsumerBase () noexcept(false)
 

Private Member Functions

void buildTrajectory (TkFittedLasBeam &beam, unsigned int &hit, vector< const GeomDetUnit *> &gd, std::vector< GlobalPoint > &globPtrack, vector< TrajectoryStateOnSurface > &tsosLas, GlobalPoint &globPref)
 
bool fitBeam (TkFittedLasBeam &beam, AlgebraicSymMatrix &covMatrix, unsigned int &hitsAtTecPlus, unsigned int &nFitParams, double &offset, double &slope, vector< GlobalPoint > &globPtrack, double &bsAngleParam, double &chi2)
 
void fitter (TkFittedLasBeam &beam, AlgebraicSymMatrix &covMatrix, unsigned int &hitsAtTecPlus, unsigned int &nFitParams, std::vector< double > &hitPhi, std::vector< double > &hitPhiError, std::vector< double > &hitZprimeError, double &zMin, double &zMax, double &bsAngleParam, double &offset, double &offsetError, double &slope, double &slopeError, double &phiAtMinusParam, double &phiAtPlusParam, double &atThetaSplitParam)
 
void getLasBeams (TkFittedLasBeam &beam, vector< TrajectoryStateOnSurface > &tsosLas)
 
void getLasHits (TkFittedLasBeam &beam, const SiStripLaserRecHit2D &hit, vector< const GeomDetUnit *> &gd, vector< GlobalPoint > &globHit, unsigned int &hitsAtTecPlus)
 
void globalTrackPoint (TkFittedLasBeam &beam, unsigned int &hit, unsigned int &hitsAtTecPlus, double &trackPhi, double &trackPhiRef, std::vector< GlobalPoint > &globHit, std::vector< GlobalPoint > &globPtrack, GlobalPoint &globPref, std::vector< double > &hitPhiError)
 
void trackPhi (TkFittedLasBeam &beam, unsigned int &hit, double &trackPhi, double &trackPhiRef, double &offset, double &slope, double &bsAngleParam, double &phiAtMinusParam, double &phiAtPlusParam, double &atThetaSplitParam, std::vector< GlobalPoint > &globHit)
 

Static Private Member Functions

static double atFunction (double *x, double *par)
 
static double tecMinusFunction (double *x, double *par)
 
static double tecPlusFunction (double *x, double *par)
 

Private Attributes

edm::ESHandle< MagneticFieldfieldHandle
 
bool fitBeamSplitters_
 
edm::Service< TFileServicefs
 
edm::ESHandle< TrackerGeometrygeometry
 
const edm::ESGetToken< TrackerGeometry, TrackerDigiGeometryRecordgeomToken_
 
TH1F * h_bsAngle
 
TH2F * h_bsAngleVsBeam
 
TH1F * h_chi2
 
TH1F * h_chi2ndof
 
TH1F * h_hitX
 
TH1F * h_hitXAt
 
TH1F * h_hitXTecMinus
 
TH1F * h_hitXTecPlus
 
TH2F * h_hitXvsZAt
 
TH2F * h_hitXvsZTecMinus
 
TH2F * h_hitXvsZTecPlus
 
TH1F * h_pull
 
TH1F * h_res
 
TH1F * h_resAt
 
TH1F * h_resTecMinus
 
TH1F * h_resTecPlus
 
TH2F * h_resVsHitAt
 
TH2F * h_resVsHitTecMinus
 
TH2F * h_resVsHitTecPlus
 
TH2F * h_resVsZAt
 
TH2F * h_resVsZTecMinus
 
TH2F * h_resVsZTecPlus
 
edm::Handle< TkLasBeamCollectionlaserBeams
 
const edm::ESGetToken< MagneticField, IdealMagneticFieldRecordmagFieldToken_
 
unsigned int nAtParameters_
 
const edm::InputTag src_
 

Static Private Attributes

static vector< double > gBarrelModuleOffset
 
static vector< double > gBarrelModuleRadius
 
static double gBeamR = 0.0
 
static double gBeamSplitterZprime = 0.0
 
static double gBeamZ0 = 0.0
 
static double gBSparam = 0.0
 
static bool gFitBeamSplitters = false
 
static unsigned int gHitsAtTecMinus = 0
 
static vector< double > gHitZprime
 
static bool gIsInnerBarrel = false
 
static float gTIBparam = 0.097614
 
static float gTOBparam = 0.034949
 

Additional Inherited Members

- Public Types inherited from edm::one::EDProducerBase
typedef EDProducerBase ModuleType
 
- Public Types inherited from edm::ProducerBase
template<typename T >
using BranchAliasSetterT = ProductRegistryHelper::BranchAliasSetterT< T >
 
using ModuleToResolverIndicies = std::unordered_multimap< std::string, std::tuple< edm::TypeID const *, const char *, edm::ProductResolverIndex > >
 
typedef ProductRegistryHelper::TypeLabelList TypeLabelList
 
- Public Types inherited from edm::EDConsumerBase
typedef ProductLabels Labels
 
- Static Public Member Functions inherited from edm::one::EDProducerBase
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &descriptions)
 
- Protected Member Functions inherited from edm::ProducerBase
template<Transition Tr = Transition::Event>
auto produces (std::string instanceName) noexcept
 declare what type of product will make and with which optional label More...
 
template<Transition B>
BranchAliasSetter produces (const TypeID &id, std::string instanceName=std::string(), bool recordProvenance=true)
 
template<BranchType B>
BranchAliasSetter produces (const TypeID &id, std::string instanceName=std::string(), bool recordProvenance=true)
 
BranchAliasSetter produces (const TypeID &id, std::string instanceName=std::string(), bool recordProvenance=true)
 
template<typename ProductType , Transition B>
BranchAliasSetterT< ProductType > produces (std::string instanceName)
 
template<class ProductType >
BranchAliasSetterT< ProductType > produces ()
 
template<typename ProductType , BranchType B>
BranchAliasSetterT< ProductType > produces (std::string instanceName)
 
template<typename ProductType , BranchType B>
BranchAliasSetterT< ProductType > produces ()
 
template<class ProductType >
BranchAliasSetterT< ProductType > produces (std::string instanceName)
 
template<typename ProductType , Transition B>
BranchAliasSetterT< ProductType > produces ()
 
template<Transition Tr = Transition::Event>
auto produces () noexcept
 
ProducesCollector producesCollector ()
 
- 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< Bconsumes (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 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

Original Authors: Gero Flucke/Kolja Kaschube Created: Wed May 6 08:43:02 CEST 2009

Id
TkLasBeamFitter.cc,v 1.12 2012/01/25 17:30:30 innocent Exp

Description: Fitting LAS beams with track model and providing TrajectoryStateOnSurface for hits.

Implementation:

Definition at line 62 of file TkLasBeamFitter.cc.

Constructor & Destructor Documentation

◆ TkLasBeamFitter()

TkLasBeamFitter::TkLasBeamFitter ( const edm::ParameterSet config)
explicit

Definition at line 210 of file TkLasBeamFitter.cc.

211  : magFieldToken_(esConsumes<edm::Transition::EndRun>()),
212  geomToken_(esConsumes<edm::Transition::EndRun>()),
213  src_(iConfig.getParameter<edm::InputTag>("src")),
214  fitBeamSplitters_(iConfig.getParameter<bool>("fitBeamSplitters")),
215  nAtParameters_(iConfig.getParameter<unsigned int>("numberOfFittedAtParameters")),
216  h_bsAngle(nullptr),
217  h_hitX(nullptr),
218  h_hitXTecPlus(nullptr),
219  h_hitXTecMinus(nullptr),
220  h_hitXAt(nullptr),
221  h_chi2(nullptr),
222  h_chi2ndof(nullptr),
223  h_pull(nullptr),
224  h_res(nullptr),
225  h_resTecPlus(nullptr),
226  h_resTecMinus(nullptr),
227  h_resAt(nullptr),
228  h_bsAngleVsBeam(nullptr),
229  h_hitXvsZTecPlus(nullptr),
230  h_hitXvsZTecMinus(nullptr),
231  h_hitXvsZAt(nullptr),
232  h_resVsZTecPlus(nullptr),
233  h_resVsZTecMinus(nullptr),
234  h_resVsZAt(nullptr),
235  h_resVsHitTecPlus(nullptr),
236  h_resVsHitTecMinus(nullptr),
237  h_resVsHitAt(nullptr) {
238  // declare the products to produce
239  this->produces<TkFittedLasBeamCollection, edm::Transition::EndRun>();
240  this->produces<TsosVectorCollection, edm::Transition::EndRun>();
241 
242  //now do what ever other initialization is needed
243 }
unsigned int nAtParameters_
const edm::ESGetToken< TrackerGeometry, TrackerDigiGeometryRecord > geomToken_
const edm::InputTag src_
const edm::ESGetToken< MagneticField, IdealMagneticFieldRecord > magFieldToken_

◆ ~TkLasBeamFitter()

TkLasBeamFitter::~TkLasBeamFitter ( )
override

Definition at line 246 of file TkLasBeamFitter.cc.

246  {
247  // do anything here that needs to be done at desctruction time
248  // (e.g. close files, deallocate resources etc.)
249 }

Member Function Documentation

◆ atFunction()

double TkLasBeamFitter::atFunction ( double *  x,
double *  par 
)
staticprivate

Definition at line 693 of file TkLasBeamFitter.cc.

References gBeamR, gBeamSplitterZprime, gBeamZ0, gBSparam, gFitBeamSplitters, gIsInnerBarrel, gTIBparam, gTOBparam, x, and z.

Referenced by fitter().

693  {
694  double z = x[0]; // 'primed'? -> yes!!!
695  // TECminus
696  if (z < -gBeamSplitterZprime - 2.0 * gBeamZ0) {
697  return par[0] + par[1] * z;
698  }
699  // BarrelMinus
700  else if (-gBeamSplitterZprime - 2.0 * gBeamZ0 < z && z < -gBeamZ0) {
701  // z value includes module offset from main beam axis
702  // TOB
703  if (!gIsInnerBarrel) {
704  return par[0] + par[1] * z + gTOBparam * (par[2] + par[4]);
705  }
706  // TIB
707  else {
708  return par[0] + par[1] * z - gTIBparam * (par[2] - par[4]);
709  }
710  }
711  // BarrelPlus
712  else if (-gBeamZ0 < z && z < gBeamSplitterZprime) {
713  // z value includes module offset from main beam axis
714  // TOB
715  if (!gIsInnerBarrel) {
716  return par[0] + par[1] * z + gTOBparam * (par[3] - par[4]);
717  }
718  // TIB
719  else {
720  return par[0] + par[1] * z - gTIBparam * (par[3] + par[4]);
721  }
722  }
723  // TECplus
724  else {
725  if (gFitBeamSplitters) {
726  // par[2] = 2*tan(BeamSplitterAngle/2.0)
727  return par[0] + par[1] * z - par[5] * (z - gBeamSplitterZprime) / gBeamR; // BS par: 5, 4, or 2
728  } else {
729  return par[0] + par[1] * z - gBSparam * (z - gBeamSplitterZprime) / gBeamR;
730  }
731  }
732 }
static float gTOBparam
static double gBeamR
static bool gFitBeamSplitters
static bool gIsInnerBarrel
static double gBeamSplitterZprime
static double gBSparam
static float gTIBparam
static double gBeamZ0

◆ buildTrajectory()

void TkLasBeamFitter::buildTrajectory ( TkFittedLasBeam beam,
unsigned int &  hit,
vector< const GeomDetUnit *> &  gd,
std::vector< GlobalPoint > &  globPtrack,
vector< TrajectoryStateOnSurface > &  tsosLas,
GlobalPoint globPref 
)
private

Definition at line 915 of file TkLasBeamFitter.cc.

References EcalCondDBWriter_cfi::beam, SurfaceSideDefinition::beforeSurface, fieldHandle, gBeamSplitterZprime, gBeamZ0, gHitZprime, HLT_2023v12_cff::magneticField, and edm::ESHandle< T >::product().

Referenced by getLasBeams().

920  {
922  GlobalVector trajectoryState;
923 
924  // TECplus
925  if (beam.isTecInternal(1)) {
926  trajectoryState = GlobalVector(globPref - globPtrack[hit]);
927  }
928  // TECminus
929  else if (beam.isTecInternal(-1)) {
930  trajectoryState = GlobalVector(globPtrack[hit] - globPref);
931  }
932  // ATs
933  else {
934  // TECminus
935  if (gHitZprime[hit] < -gBeamSplitterZprime - 2.0 * gBeamZ0) {
936  trajectoryState = GlobalVector(globPtrack[hit] - globPref);
937  }
938  // TECplus
939  else if (gHitZprime[hit] > gBeamSplitterZprime) {
940  trajectoryState = GlobalVector(globPref - globPtrack[hit]);
941  }
942  // Barrel
943  else {
944  trajectoryState = GlobalVector(globPtrack[hit] - globPref);
945  }
946  }
947  // cout << "trajectory: " << trajectoryState << endl;
948  const FreeTrajectoryState ftsLas = FreeTrajectoryState(globPtrack[hit], trajectoryState, 0, magneticField);
949  tsosLas.push_back(TrajectoryStateOnSurface(ftsLas, gd[hit]->surface(), SurfaceSideDefinition::beforeSurface));
950 }
static vector< double > gHitZprime
edm::ESHandle< MagneticField > fieldHandle
T const * product() const
Definition: ESHandle.h:86
static double gBeamSplitterZprime
Global3DVector GlobalVector
Definition: GlobalVector.h:10
static double gBeamZ0

◆ endRunProduce()

void TkLasBeamFitter::endRunProduce ( edm::Run run,
const edm::EventSetup setup 
)
override

Definition at line 263 of file TkLasBeamFitter.cc.

References EcalCondDBWriter_cfi::beam, gather_cfg::cout, fieldHandle, fitBeamSplitters_, fs, gBSparam, geomToken_, getLasBeams(), gFitBeamSplitters, h_bsAngle, h_bsAngleVsBeam, h_chi2, h_chi2ndof, h_hitX, h_hitXAt, h_hitXTecMinus, h_hitXTecPlus, h_hitXvsZAt, h_hitXvsZTecMinus, h_hitXvsZTecPlus, h_pull, h_res, h_resAt, h_resTecMinus, h_resTecPlus, h_resVsHitAt, h_resVsHitTecMinus, h_resVsHitTecPlus, h_resVsZAt, h_resVsZTecMinus, h_resVsZTecPlus, laserBeams, magFieldToken_, TFileService::make(), eostools::move(), writedatasetfile::run, and singleTopDQM_cfi::setup.

263  {
264  // }
265  // // FIXME!
266  // // Indeed, that should be in endRun(..) - as soon as AlignmentProducer can call
267  // // the algorithm's endRun correctly!
268  //
269  //
270  // void TkLasBeamFitter::beginRun(edm::Run &run, const edm::EventSetup &setup)
271  // {
272 
273  // book histograms
274  h_hitX = fs->make<TH1F>("hitX", "local x of LAS hits;local x [cm];N", 100, -0.5, 0.5);
275  h_hitXTecPlus = fs->make<TH1F>("hitXTecPlus", "local x of LAS hits in TECplus;local x [cm];N", 100, -0.5, 0.5);
276  h_hitXTecMinus = fs->make<TH1F>("hitXTecMinus", "local x of LAS hits in TECminus;local x [cm];N", 100, -0.5, 0.5);
277  h_hitXAt = fs->make<TH1F>("hitXAt", "local x of LAS hits in ATs;local x [cm];N", 100, -2.5, 2.5);
279  fs->make<TH2F>("hitXvsZTecPlus", "local x vs z in TECplus;z [cm];local x [cm]", 80, 120, 280, 100, -0.5, 0.5);
281  fs->make<TH2F>("hitXvsZTecMinus", "local x vs z in TECMinus;z [cm];local x [cm]", 80, -280, -120, 100, -0.5, 0.5);
282  h_hitXvsZAt = fs->make<TH2F>("hitXvsZAt", "local x vs z in ATs;z [cm];local x [cm]", 200, -200, 200, 100, -0.5, 0.5);
283  h_chi2 = fs->make<TH1F>("chi2", "#chi^{2};#chi^{2};N", 100, 0, 2000);
284  h_chi2ndof = fs->make<TH1F>("chi2ndof", "#chi^{2} per degree of freedom;#chi^{2}/N_{dof};N", 100, 0, 300);
285  h_pull = fs->make<TH1F>("pull", "pulls of #phi residuals;pull;N", 50, -10, 10);
286  h_res = fs->make<TH1F>("res", "#phi residuals;#phi_{track} - #phi_{hit} [rad];N", 60, -0.0015, 0.0015);
287  h_resTecPlus =
288  fs->make<TH1F>("resTecPlus", "#phi residuals in TECplus;#phi_{track} - #phi_{hit} [rad];N", 30, -0.0015, 0.0015);
289  h_resTecMinus = fs->make<TH1F>(
290  "resTecMinus", "#phi residuals in TECminus;#phi_{track} - #phi_{hit} [rad];N", 60, -0.0015, 0.0015);
291  h_resAt = fs->make<TH1F>("resAt", "#phi residuals in ATs;#phi_{track} - #phi_{hit} [rad];N", 30, -0.0015, 0.0015);
292  h_resVsZTecPlus = fs->make<TH2F>("resVsZTecPlus",
293  "phi residuals vs. z in TECplus;z [cm];#phi_{track} - #phi_{hit} [rad]",
294  80,
295  120,
296  280,
297  100,
298  -0.0015,
299  0.0015);
300  h_resVsZTecMinus = fs->make<TH2F>("resVsZTecMinus",
301  "phi residuals vs. z in TECminus;z [cm];#phi_{track} - #phi_{hit} [rad]",
302  80,
303  -280,
304  -120,
305  100,
306  -0.0015,
307  0.0015);
308  h_resVsZAt = fs->make<TH2F>(
309  "resVsZAt", "#phi residuals vs. z in ATs;N;#phi_{track} - #phi_{hit} [rad]", 200, -200, 200, 100, -0.0015, 0.0015);
310  h_resVsHitTecPlus = fs->make<TH2F>("resVsHitTecPlus",
311  "#phi residuals vs. hits in TECplus;hit no.;#phi_{track} - #phi_{hit} [rad]",
312  144,
313  0,
314  144,
315  100,
316  -0.0015,
317  0.0015);
318  h_resVsHitTecMinus = fs->make<TH2F>("resVsHitTecMinus",
319  "#phi residuals vs. hits in TECminus;hit no.;#phi_{track} - #phi_{hit} [rad]",
320  144,
321  0,
322  144,
323  100,
324  -0.0015,
325  0.0015);
326  h_resVsHitAt = fs->make<TH2F>("resVsHitAt",
327  "#phi residuals vs. hits in ATs;hit no.;#phi_{track} - #phi_{hit} [rad]",
328  176,
329  0,
330  176,
331  100,
332  -0.0015,
333  0.0015);
334  h_bsAngle = fs->make<TH1F>("bsAngle", "fitted beam splitter angle;BS angle [rad];N", 40, -0.004, 0.004);
335  h_bsAngleVsBeam = fs->make<TH2F>(
336  "bsAngleVsBeam", "fitted beam splitter angle per beam;Beam no.;BS angle [rad]", 40, 0, 300, 100, -0.004, 0.004);
337 
338  // Create output collections - they are parallel.
339  // (edm::Ref etc. and thus edm::AssociationVector are not supported for edm::Run...)
340  auto fittedBeams = std::make_unique<TkFittedLasBeamCollection>();
341  // One std::vector<TSOS> for each TkFittedLasBeam:
342  auto tsosesVec = std::make_unique<TsosVectorCollection>();
343 
344  // get TkLasBeams, Tracker geometry, magnetic field
345  run.getByLabel("LaserAlignment", "tkLaserBeams", laserBeams);
346  geometry = setup.getHandle(geomToken_);
347  fieldHandle = setup.getHandle(magFieldToken_);
348 
349  // hack for fixed BSparams (ugly!)
350  // double bsParams[34] = {-0.000266,-0.000956,-0.001205,-0.000018,-0.000759,0.002554,
351  // 0.000465,0.000975,0.001006,0.002027,-0.001263,-0.000763,
352  // -0.001702,0.000906,-0.002120,0.001594,0.000661,-0.000457,
353  // -0.000447,0.000347,-0.002266,-0.000446,0.000659,0.000018,
354  // -0.001630,-0.000324,
355  // // ATs
356  // -999.,-0.001709,-0.002091,-999.,
357  // -0.001640,-999.,-0.002444,-0.002345};
358 
359  double bsParams[40] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
360  0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
361 
362  // beam counter
363  unsigned int beamNo(0);
364  // fit BS? If false, values from bsParams are taken
366  if (fitBeamSplitters_)
367  cout << "Fitting BS!" << endl;
368  else
369  cout << "BS fixed, not fitted!" << endl;
370 
371  // loop over LAS beams
372  for (TkLasBeamCollection::const_iterator iBeam = laserBeams->begin(), iEnd = laserBeams->end(); iBeam != iEnd;
373  ++iBeam) {
374  TkFittedLasBeam beam(*iBeam);
375  vector<TrajectoryStateOnSurface> tsosLas;
376 
377  // set BS param for fit
378  gBSparam = bsParams[beamNo];
379 
380  // call main function; all other functions are called inside getLasBeams(..)
381  this->getLasBeams(beam, tsosLas);
382 
383  // fill output products
384  fittedBeams->push_back(beam);
385  tsosesVec->push_back(tsosLas);
386 
387  // if(!this->fitBeam(fittedBeams->back(), tsosesVec->back())){
388  // edm::LogError("BadFit")
389  // << "Problems fitting TkLasBeam, id " << fittedBeams->back().getBeamId() << ".";
390  // fittedBeams->pop_back(); // remove last entry added just before
391  // tsosesVec->pop_back(); // dito
392  // }
393 
394  beamNo++;
395  }
396 
397  // finally, put fitted beams and TSOS vectors into run
398  run.put(std::move(fittedBeams));
399  run.put(std::move(tsosesVec));
400 }
static bool gFitBeamSplitters
edm::ESHandle< MagneticField > fieldHandle
edm::Handle< TkLasBeamCollection > laserBeams
const edm::ESGetToken< TrackerGeometry, TrackerDigiGeometryRecord > geomToken_
static double gBSparam
const edm::ESGetToken< MagneticField, IdealMagneticFieldRecord > magFieldToken_
T * make(const Args &...args) const
make new ROOT object
Definition: TFileService.h:64
void getLasBeams(TkFittedLasBeam &beam, vector< TrajectoryStateOnSurface > &tsosLas)
edm::Service< TFileService > fs
def move(src, dest)
Definition: eostools.py:511

◆ fitBeam()

bool TkLasBeamFitter::fitBeam ( TkFittedLasBeam beam,
AlgebraicSymMatrix covMatrix,
unsigned int &  hitsAtTecPlus,
unsigned int &  nFitParams,
double &  offset,
double &  slope,
vector< GlobalPoint > &  globPtrack,
double &  bsAngleParam,
double &  chi2 
)
private

Definition at line 953 of file TkLasBeamFitter.cc.

References EcalCondDBWriter_cfi::beam, nano_mu_local_reco_cff::chi2, fitBeamSplitters_, gBarrelModuleOffset, gBeamSplitterZprime, gBeamZ0, gHitsAtTecMinus, gHitZprime, hltrates_dqm_sourceclient-live_cfg::offset, submitPVValidationJobs::params, and slope.

Referenced by getLasBeams().

961  {
962  // set beam parameters for beam output
963  unsigned int paramType(0);
964  if (!fitBeamSplitters_)
965  paramType = 1;
966  if (beam.isAlignmentTube() && hitsAtTecPlus == 0)
967  paramType = 0;
968  // const unsigned int nPedeParams = nFitParams + paramType;
969 
970  // test without BS params
971  const unsigned int nPedeParams(nFitParams);
972  // cout << "number of Pede parameters: " << nPedeParams << endl;
973 
974  std::vector<TkFittedLasBeam::Scalar> params(nPedeParams);
975  params[0] = offset;
976  params[1] = slope;
977  // no BS parameter for AT beams without TECplus hits
978  // if(beam.isTecInternal() || hitsAtTecPlus > 0) params[2] = bsAngleParam;
979 
980  AlgebraicMatrix derivatives(gHitZprime.size(), nPedeParams);
981  // fill derivatives matrix with local track derivatives
982  for (unsigned int hit = 0; hit < gHitZprime.size(); ++hit) {
983  // d(delta phi)/d(offset) is identical for every hit
984  derivatives[hit][0] = 1.0;
985 
986  // d(delta phi)/d(slope) and d(delta phi)/d(bsAngleParam) depend on parametrizations
987  // TECplus
988  if (beam.isTecInternal(1)) {
989  derivatives[hit][1] = globPtrack[hit].z();
990  // if(gHitZprime[hit] < gBeamSplitterZprime){
991  // derivatives[hit][2] = 0.0;
992  // }
993  // else{
994  // derivatives[hit][2] = - (globPtrack[hit].z() - gBeamSplitterZprime) / gBeamR;
995  // }
996  }
997  // TECminus
998  else if (beam.isTecInternal(-1)) {
999  derivatives[hit][1] = globPtrack[hit].z();
1000  // if(gHitZprime[hit] > gBeamSplitterZprime){
1001  // derivatives[hit][2] = 0.0;
1002  // }
1003  // else{
1004  // derivatives[hit][2] = (globPtrack[hit].z() - gBeamSplitterZprime) / gBeamR;
1005  // }
1006  }
1007  // ATs
1008  else {
1009  // TECminus
1010  if (gHitZprime[hit] < -gBeamSplitterZprime - 2.0 * gBeamZ0) {
1011  derivatives[hit][1] = globPtrack[hit].z();
1012  // if(hitsAtTecPlus > 0){
1013  // derivatives[hit][2] = 0.0;
1014  // }
1015  }
1016  // TECplus
1017  else if (gHitZprime[hit] > gBeamSplitterZprime) {
1018  derivatives[hit][1] = globPtrack[hit].z();
1019  // if(hitsAtTecPlus > 0){
1020  // derivatives[hit][2] = - (globPtrack[hit].z() - gBeamSplitterZprime) / gBeamR;
1021  // }
1022  }
1023  // Barrel
1024  else {
1025  derivatives[hit][1] = globPtrack[hit].z() - gBarrelModuleOffset[hit - gHitsAtTecMinus];
1026  // if(hitsAtTecPlus > 0){
1027  // derivatives[hit][2] = 0.0;
1028  // }
1029  }
1030  }
1031  }
1032 
1033  unsigned int firstFixedParam(covMatrix.num_col()); // FIXME! --> no, is fine!!!
1034  // unsigned int firstFixedParam = nPedeParams - 1;
1035  // if(beam.isAlignmentTube() && hitsAtTecPlus == 0) firstFixedParam = nPedeParams;
1036  // cout << "first fixed parameter: " << firstFixedParam << endl;
1037  // set fit results
1038  beam.setParameters(paramType, params, covMatrix, derivatives, firstFixedParam, chi2);
1039 
1040  return true; // return false in case of problems
1041 }
static vector< double > gHitZprime
static const double slope[3]
CLHEP::HepMatrix AlgebraicMatrix
static double gBeamSplitterZprime
static unsigned int gHitsAtTecMinus
static vector< double > gBarrelModuleOffset
static double gBeamZ0

◆ fitter()

void TkLasBeamFitter::fitter ( TkFittedLasBeam beam,
AlgebraicSymMatrix covMatrix,
unsigned int &  hitsAtTecPlus,
unsigned int &  nFitParams,
std::vector< double > &  hitPhi,
std::vector< double > &  hitPhiError,
std::vector< double > &  hitZprimeError,
double &  zMin,
double &  zMax,
double &  bsAngleParam,
double &  offset,
double &  offsetError,
double &  slope,
double &  slopeError,
double &  phiAtMinusParam,
double &  phiAtPlusParam,
double &  atThetaSplitParam 
)
private

Definition at line 735 of file TkLasBeamFitter.cc.

References atFunction(), EcalCondDBWriter_cfi::beam, cuy::col, fitBeamSplitters_, gBSparam, gHitZprime, hltrates_dqm_sourceclient-live_cfg::offset, slope, tecMinusFunction(), tecPlusFunction(), photonAnalyzer_cfi::zMax, and photonAnalyzer_cfi::zMin.

Referenced by getLasBeams().

751  {
752  TGraphErrors *lasData =
753  new TGraphErrors(gHitZprime.size(), &(gHitZprime[0]), &(hitPhi[0]), &(hitZprimeError[0]), &(hitPhiError[0]));
754 
755  // do fit (R = entire range)
756  if (beam.isTecInternal(1)) {
757  TF1 tecPlus("tecPlus", tecPlusFunction, zMin, zMax, nFitParams);
758  tecPlus.SetParameter(1, 0); // slope
759  tecPlus.SetParameter(nFitParams - 1, 0); // BS
760  lasData->Fit(&tecPlus, "R"); // "R", "RV" or "RQ"
761  } else if (beam.isTecInternal(-1)) {
762  TF1 tecMinus("tecMinus", tecMinusFunction, zMin, zMax, nFitParams);
763  tecMinus.SetParameter(1, 0); // slope
764  tecMinus.SetParameter(nFitParams - 1, 0); // BS
765  lasData->Fit(&tecMinus, "R");
766  } else {
767  TF1 at("at", atFunction, zMin, zMax, nFitParams);
768  at.SetParameter(1, 0); // slope
769  at.SetParameter(nFitParams - 1, 0); // BS
770  lasData->Fit(&at, "R");
771  }
772 
773  // get values and errors for offset and slope
774  gMinuit->GetParameter(0, offset, offsetError);
775  gMinuit->GetParameter(1, slope, slopeError);
776 
777  // additional AT parameters
778  // define param errors that are not used later
779  double bsAngleParamError(0.), phiAtMinusParamError(0.), phiAtPlusParamError(0.), atThetaSplitParamError(0.);
780 
781  if (beam.isAlignmentTube()) {
782  gMinuit->GetParameter(2, phiAtMinusParam, phiAtMinusParamError);
783  gMinuit->GetParameter(3, phiAtPlusParam, phiAtPlusParamError);
784  gMinuit->GetParameter(4, atThetaSplitParam, atThetaSplitParamError);
785  }
786  // get Beam Splitter parameters
787  if (fitBeamSplitters_) {
788  if (beam.isAlignmentTube() && hitsAtTecPlus == 0) {
789  bsAngleParam = gBSparam;
790  } else {
791  gMinuit->GetParameter(nFitParams - 1, bsAngleParam, bsAngleParamError);
792  }
793  } else {
794  bsAngleParam = gBSparam;
795  }
796 
797  // fill covariance matrix
798  vector<double> vec(covMatrix.num_col() * covMatrix.num_col());
799  gMinuit->mnemat(&vec[0], covMatrix.num_col());
800  for (int col = 0; col < covMatrix.num_col(); col++) {
801  for (int row = 0; row < covMatrix.num_col(); row++) {
802  covMatrix[col][row] = vec[row + covMatrix.num_col() * col];
803  }
804  }
805  // compute correlation between parameters
806  // double corr01 = covMatrix[1][0]/(offsetError*slopeError);
807 
808  delete lasData;
809 }
static vector< double > gHitZprime
static double tecMinusFunction(double *x, double *par)
static double tecPlusFunction(double *x, double *par)
static const double slope[3]
static double gBSparam
col
Definition: cuy.py:1009
static double atFunction(double *x, double *par)

◆ getLasBeams()

void TkLasBeamFitter::getLasBeams ( TkFittedLasBeam beam,
vector< TrajectoryStateOnSurface > &  tsosLas 
)
private

Fit 'beam' using info from its base class TkLasBeam and set its parameters. Also fill 'tsoses' with TSOS for each LAS hit.

Definition at line 405 of file TkLasBeamFitter.cc.

References funct::abs(), EcalCondDBWriter_cfi::beam, buildTrajectory(), nano_mu_local_reco_cff::chi2, gather_cfg::cout, fitBeam(), fitBeamSplitters_, fitter(), gBarrelModuleOffset, gBarrelModuleRadius, gBeamR, gBeamSplitterZprime, gBeamZ0, getLasHits(), gHitsAtTecMinus, gHitZprime, gIsInnerBarrel, globalTrackPoint(), h_bsAngle, h_bsAngleVsBeam, h_chi2, h_chi2ndof, h_hitX, h_hitXAt, h_hitXTecMinus, h_hitXTecPlus, h_hitXvsZAt, h_hitXvsZTecMinus, h_hitXvsZTecPlus, h_pull, h_res, h_resAt, h_resTecMinus, h_resTecPlus, h_resVsHitAt, h_resVsHitTecMinus, h_resVsHitTecPlus, h_resVsZAt, h_resVsZTecMinus, h_resVsZTecPlus, nAtParameters_, hltrates_dqm_sourceclient-live_cfg::offset, perp(), phi, slope, trackPhi(), hit::x, z, photonAnalyzer_cfi::zMax, and photonAnalyzer_cfi::zMin.

Referenced by endRunProduce().

405  {
406  cout << "---------------------------------------" << endl;
407  cout << "beam id: " << beam.getBeamId() // << " isTec: " << (beam.isTecInternal() ? "Y" : "N")
408  << " isTec+: " << (beam.isTecInternal(1) ? "Y" : "N") << " isTec-: " << (beam.isTecInternal(-1) ? "Y" : "N")
409  << " isAt: " << (beam.isAlignmentTube() ? "Y" : "N") << " isR6: " << (beam.isRing6() ? "Y" : "N") << endl;
410 
411  // reset static variables
412  gHitsAtTecMinus = 0;
413  gHitZprime.clear();
414  gBarrelModuleRadius.clear();
415  gBarrelModuleOffset.clear();
416 
417  // set right beam radius
418  gBeamR = beam.isRing6() ? 84.0 : 56.4;
419 
420  vector<const GeomDetUnit *> gd;
421  vector<GlobalPoint> globHit;
422  unsigned int hitsAtTecPlus(0);
423  double sumZ(0.);
424 
425  // loop over hits
426  for (TkLasBeam::const_iterator iHit = beam.begin(); iHit < beam.end(); ++iHit) {
427  // iHit is a SiStripLaserRecHit2D
428 
429  const SiStripLaserRecHit2D &hit(*iHit);
430 
431  this->getLasHits(beam, hit, gd, globHit, hitsAtTecPlus);
432  sumZ += globHit.back().z();
433 
434  // fill histos
435  h_hitX->Fill(hit.localPosition().x());
436  // TECplus
437  if (beam.isTecInternal(1)) {
438  h_hitXTecPlus->Fill(hit.localPosition().x());
439  h_hitXvsZTecPlus->Fill(globHit.back().z(), hit.localPosition().x());
440  }
441  // TECminus
442  else if (beam.isTecInternal(-1)) {
443  h_hitXTecMinus->Fill(hit.localPosition().x());
444  h_hitXvsZTecMinus->Fill(globHit.back().z(), hit.localPosition().x());
445  }
446  // ATs
447  else {
448  h_hitXAt->Fill(hit.localPosition().x());
449  h_hitXvsZAt->Fill(globHit.back().z(), hit.localPosition().x());
450  }
451  }
452 
453  gBeamZ0 = sumZ / globHit.size();
454  double zMin(0.), zMax(0.);
455  // TECplus
456  if (beam.isTecInternal(1)) {
457  gBeamSplitterZprime = 205.75 - gBeamZ0;
458  zMin = 120.0 - gBeamZ0;
459  zMax = 280.0 - gBeamZ0;
460  }
461  // TECminus
462  else if (beam.isTecInternal(-1)) {
463  gBeamSplitterZprime = -205.75 - gBeamZ0;
464  zMin = -280.0 - gBeamZ0;
465  zMax = -120.0 - gBeamZ0;
466  }
467  // AT
468  else {
469  gBeamSplitterZprime = 112.3 - gBeamZ0;
470  zMin = -200.0 - gBeamZ0;
471  zMax = 200.0 - gBeamZ0;
472  }
473 
474  // fill vectors for fitted quantities
475  vector<double> hitPhi, hitPhiError, hitZprimeError;
476 
477  for (unsigned int hit = 0; hit < globHit.size(); ++hit) {
478  hitPhi.push_back(static_cast<double>(globHit[hit].phi()));
479  // localPositionError[hit] or assume 0.003, 0.006
480  hitPhiError.push_back(0.003 / globHit[hit].perp());
481  // no errors on z, fill with zeros
482  hitZprimeError.push_back(0.0);
483  // barrel-specific values
484  if (beam.isAlignmentTube() && abs(globHit[hit].z()) < 112.3) {
485  gBarrelModuleRadius.push_back(globHit[hit].perp());
486  gBarrelModuleOffset.push_back(gBarrelModuleRadius.back() - gBeamR);
487  // TIB/TOB flag
488  if (gBarrelModuleOffset.back() < 0.0) {
489  gIsInnerBarrel = true;
490  } else {
491  gIsInnerBarrel = false;
492  }
493  gHitZprime.push_back(globHit[hit].z() - gBeamZ0 - abs(gBarrelModuleOffset.back()));
494  }
495  // non-barrel z'
496  else {
497  gHitZprime.push_back(globHit[hit].z() - gBeamZ0);
498  }
499  }
500 
501  // number of fit parameters, 3 for TECs (always!); 3, 5, or 6 for ATs
502  unsigned int tecParams(3), atParams(0);
503  if (nAtParameters_ == 3)
504  atParams = 3;
505  else if (nAtParameters_ == 5)
506  atParams = 5;
507  else
508  atParams = 6; // <-- default value, recommended
509  unsigned int nFitParams(0);
510  if (!fitBeamSplitters_ || (hitsAtTecPlus == 0 && beam.isAlignmentTube())) {
511  tecParams = tecParams - 1;
512  atParams = atParams - 1;
513  }
514  if (beam.isTecInternal()) {
515  nFitParams = tecParams;
516  } else {
517  nFitParams = atParams;
518  }
519 
520  // fit parameter definitions
521  double offset(0.), offsetError(0.), slope(0.), slopeError(0.), bsAngleParam(0.), phiAtMinusParam(0.),
522  phiAtPlusParam(0.), atThetaSplitParam(0.);
523  AlgebraicSymMatrix covMatrix;
524  if (!fitBeamSplitters_ || (beam.isAlignmentTube() && hitsAtTecPlus == 0)) {
525  covMatrix = AlgebraicSymMatrix(nFitParams, 1);
526  } else {
527  covMatrix = AlgebraicSymMatrix(nFitParams - 1, 1);
528  }
529 
530  this->fitter(beam,
531  covMatrix,
532  hitsAtTecPlus,
533  nFitParams,
534  hitPhi,
535  hitPhiError,
536  hitZprimeError,
537  zMin,
538  zMax,
539  bsAngleParam,
540  offset,
541  offsetError,
542  slope,
543  slopeError,
544  phiAtMinusParam,
545  phiAtPlusParam,
546  atThetaSplitParam);
547 
548  vector<GlobalPoint> globPtrack;
549  GlobalPoint globPref;
550  double chi2(0.);
551 
552  for (unsigned int hit = 0; hit < gHitZprime.size(); ++hit) {
553  // additional phi value (trackPhiRef) for trajectory calculation
554  double trackPhi(0.), trackPhiRef(0.);
555 
556  this->trackPhi(beam,
557  hit,
558  trackPhi,
559  trackPhiRef,
560  offset,
561  slope,
562  bsAngleParam,
563  phiAtMinusParam,
564  phiAtPlusParam,
565  atThetaSplitParam,
566  globHit);
567 
568  cout << "track phi = " << trackPhi << ", hit phi = " << hitPhi[hit] << ", zPrime = " << gHitZprime[hit]
569  << " r = " << globHit[hit].perp() << endl;
570 
571  this->globalTrackPoint(beam, hit, hitsAtTecPlus, trackPhi, trackPhiRef, globHit, globPtrack, globPref, hitPhiError);
572 
573  // calculate residuals = pred - hit (in global phi)
574  const double phiResidual = globPtrack[hit].phi() - globHit[hit].phi();
575  // pull calculation (FIX!!!)
576  const double phiResidualPull = phiResidual / hitPhiError[hit];
577  // sqrt(hitPhiError[hit]*hitPhiError[hit] +
578  // (offsetError*offsetError + globPtrack[hit].z()*globPtrack[hit].z() * slopeError*slopeError));
579 
580  // calculate chi2
581  chi2 += phiResidual * phiResidual / (hitPhiError[hit] * hitPhiError[hit]);
582 
583  // fill histos
584  h_res->Fill(phiResidual);
585  // TECplus
586  if (beam.isTecInternal(1)) {
587  h_pull->Fill(phiResidualPull);
588  h_resTecPlus->Fill(phiResidual);
589  h_resVsZTecPlus->Fill(globPtrack[hit].z(), phiResidual);
590  // Ring 6
591  if (beam.isRing6()) {
592  h_resVsHitTecPlus->Fill(hit + (beam.getBeamId() - 1) / 10 * 9 + 72, phiResidual);
593  }
594  // Ring 4
595  else {
596  h_resVsHitTecPlus->Fill(hit + beam.getBeamId() / 10 * 9, phiResidual);
597  }
598  }
599  // TECminus
600  else if (beam.isTecInternal(-1)) {
601  h_pull->Fill(phiResidualPull);
602  h_resTecMinus->Fill(phiResidual);
603  h_resVsZTecMinus->Fill(globPtrack[hit].z(), phiResidual);
604  // Ring 6
605  if (beam.isRing6()) {
606  h_resVsHitTecMinus->Fill(hit + (beam.getBeamId() - 101) / 10 * 9 + 72, phiResidual);
607  }
608  // Ring 4
609  else {
610  h_resVsHitTecMinus->Fill(hit + (beam.getBeamId() - 100) / 10 * 9, phiResidual);
611  }
612  }
613  // ATs
614  else {
615  h_pull->Fill(phiResidualPull);
616  h_resAt->Fill(phiResidual);
617  h_resVsZAt->Fill(globPtrack[hit].z(), phiResidual);
618  h_resVsHitAt->Fill(hit + (beam.getBeamId() - 200) / 10 * 22, phiResidual);
619  }
620 
621  this->buildTrajectory(beam, hit, gd, globPtrack, tsosLas, globPref);
622  }
623 
624  cout << "chi^2 = " << chi2 << ", chi^2/ndof = " << chi2 / (gHitZprime.size() - nFitParams) << endl;
625  this->fitBeam(beam, covMatrix, hitsAtTecPlus, nFitParams, offset, slope, globPtrack, bsAngleParam, chi2);
626 
627  cout << "bsAngleParam = " << bsAngleParam << endl;
628 
629  // fill histos
630  // include slope, offset, covariance plots here
631  h_chi2->Fill(chi2);
632  h_chi2ndof->Fill(chi2 / (gHitZprime.size() - nFitParams));
633  if (bsAngleParam != 0.0) {
634  h_bsAngle->Fill(2.0 * atan(0.5 * bsAngleParam));
635  h_bsAngleVsBeam->Fill(beam.getBeamId(), 2.0 * atan(0.5 * bsAngleParam));
636  }
637 }
static vector< double > gHitZprime
static double gBeamR
static const double slope[3]
static vector< double > gBarrelModuleRadius
bool fitBeam(TkFittedLasBeam &beam, AlgebraicSymMatrix &covMatrix, unsigned int &hitsAtTecPlus, unsigned int &nFitParams, double &offset, double &slope, vector< GlobalPoint > &globPtrack, double &bsAngleParam, double &chi2)
void getLasHits(TkFittedLasBeam &beam, const SiStripLaserRecHit2D &hit, vector< const GeomDetUnit *> &gd, vector< GlobalPoint > &globHit, unsigned int &hitsAtTecPlus)
static bool gIsInnerBarrel
void buildTrajectory(TkFittedLasBeam &beam, unsigned int &hit, vector< const GeomDetUnit *> &gd, std::vector< GlobalPoint > &globPtrack, vector< TrajectoryStateOnSurface > &tsosLas, GlobalPoint &globPref)
unsigned int nAtParameters_
void globalTrackPoint(TkFittedLasBeam &beam, unsigned int &hit, unsigned int &hitsAtTecPlus, double &trackPhi, double &trackPhiRef, std::vector< GlobalPoint > &globHit, std::vector< GlobalPoint > &globPtrack, GlobalPoint &globPref, std::vector< double > &hitPhiError)
static double gBeamSplitterZprime
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
T perp() const
Magnitude of transverse component.
static unsigned int gHitsAtTecMinus
void fitter(TkFittedLasBeam &beam, AlgebraicSymMatrix &covMatrix, unsigned int &hitsAtTecPlus, unsigned int &nFitParams, std::vector< double > &hitPhi, std::vector< double > &hitPhiError, std::vector< double > &hitZprimeError, double &zMin, double &zMax, double &bsAngleParam, double &offset, double &offsetError, double &slope, double &slopeError, double &phiAtMinusParam, double &phiAtPlusParam, double &atThetaSplitParam)
static vector< double > gBarrelModuleOffset
void trackPhi(TkFittedLasBeam &beam, unsigned int &hit, double &trackPhi, double &trackPhiRef, double &offset, double &slope, double &bsAngleParam, double &phiAtMinusParam, double &phiAtPlusParam, double &atThetaSplitParam, std::vector< GlobalPoint > &globHit)
CLHEP::HepSymMatrix AlgebraicSymMatrix
std::vector< SiStripLaserRecHit2D >::const_iterator const_iterator
Definition: TkLasBeam.h:14
static double gBeamZ0

◆ getLasHits()

void TkLasBeamFitter::getLasHits ( TkFittedLasBeam beam,
const SiStripLaserRecHit2D hit,
vector< const GeomDetUnit *> &  gd,
vector< GlobalPoint > &  globHit,
unsigned int &  hitsAtTecPlus 
)
private

Definition at line 640 of file TkLasBeamFitter.cc.

References funct::abs(), EcalCondDBWriter_cfi::beam, gHitsAtTecMinus, and PV3DBase< T, PVType, FrameType >::z().

Referenced by getLasBeams().

644  {
645  // get global position of LAS hits
646  gd.push_back(geometry->idToDetUnit(hit.getDetId()));
647  GlobalPoint globPtemp(gd.back()->toGlobal(hit.localPosition()));
648 
649  // testing: globPtemp should be right
650  globHit.push_back(globPtemp);
651 
652  if (beam.isAlignmentTube()) {
653  if (abs(globPtemp.z()) > 112.3) {
654  if (globPtemp.z() < 112.3)
655  gHitsAtTecMinus++;
656  else
657  hitsAtTecPlus++;
658  }
659  }
660 }
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
static unsigned int gHitsAtTecMinus

◆ globalTrackPoint()

void TkLasBeamFitter::globalTrackPoint ( TkFittedLasBeam beam,
unsigned int &  hit,
unsigned int &  hitsAtTecPlus,
double &  trackPhi,
double &  trackPhiRef,
std::vector< GlobalPoint > &  globHit,
std::vector< GlobalPoint > &  globPtrack,
GlobalPoint globPref,
std::vector< double > &  hitPhiError 
)
private

Definition at line 880 of file TkLasBeamFitter.cc.

References EcalCondDBWriter_cfi::beam, gBeamR, gHitsAtTecMinus, gHitZprime, perp(), trackPhi(), and z.

Referenced by getLasBeams().

888  {
889  // TECs
890  if (beam.isTecInternal(0)) {
891  globPtrack.push_back(GlobalPoint(GlobalPoint::Cylindrical(gBeamR, trackPhi, globHit[hit].z())));
892  globPref = GlobalPoint(GlobalPoint::Cylindrical(gBeamR, trackPhiRef, globHit[hit].z() + 1.0));
893  }
894  // ATs
895  else {
896  // TECminus
897  if (hit < gHitsAtTecMinus) { // gHitZprime[hit] < -gBeamSplitterZprime - 2.0*gBeamZ0
898  globPtrack.push_back(GlobalPoint(GlobalPoint::Cylindrical(gBeamR, trackPhi, globHit[hit].z())));
899  globPref = GlobalPoint(GlobalPoint::Cylindrical(gBeamR, trackPhiRef, globHit[hit].z() + 1.0));
900  }
901  // TECplus
902  else if (hit > gHitZprime.size() - hitsAtTecPlus - 1) { // gHitZprime[hit] > gBeamSplitterZprime
903  globPtrack.push_back(GlobalPoint(GlobalPoint::Cylindrical(gBeamR, trackPhi, globHit[hit].z())));
904  globPref = GlobalPoint(GlobalPoint::Cylindrical(gBeamR, trackPhiRef, globHit[hit].z() + 1.0));
905  }
906  // Barrel
907  else {
908  globPtrack.push_back(GlobalPoint(GlobalPoint::Cylindrical(globHit[hit].perp(), trackPhi, globHit[hit].z())));
909  globPref = GlobalPoint(GlobalPoint::Cylindrical(gBeamR, trackPhiRef, globHit[hit].z()));
910  }
911  }
912 }
static vector< double > gHitZprime
static double gBeamR
Global3DPoint GlobalPoint
Definition: GlobalPoint.h:10
T perp() const
Magnitude of transverse component.
static unsigned int gHitsAtTecMinus
void trackPhi(TkFittedLasBeam &beam, unsigned int &hit, double &trackPhi, double &trackPhiRef, double &offset, double &slope, double &bsAngleParam, double &phiAtMinusParam, double &phiAtPlusParam, double &atThetaSplitParam, std::vector< GlobalPoint > &globHit)

◆ produce()

void TkLasBeamFitter::produce ( edm::Event event,
const edm::EventSetup setup 
)
overridevirtual

Implements edm::one::EDProducerBase.

Definition at line 257 of file TkLasBeamFitter.cc.

257  {
258  // Nothing per event!
259 }

◆ tecMinusFunction()

double TkLasBeamFitter::tecMinusFunction ( double *  x,
double *  par 
)
staticprivate

Definition at line 678 of file TkLasBeamFitter.cc.

References gBeamR, gBeamSplitterZprime, gBSparam, gFitBeamSplitters, x, and z.

Referenced by fitter().

678  {
679  double z = x[0]; // 'primed'? -> yes!!!
680 
681  if (z > gBeamSplitterZprime) {
682  return par[0] + par[1] * z;
683  } else {
684  if (gFitBeamSplitters) {
685  // par[2] = 2*tan(BeamSplitterAngle/2.0)
686  return par[0] + par[1] * z + par[2] * (z - gBeamSplitterZprime) / gBeamR;
687  } else {
688  return par[0] + par[1] * z + gBSparam * (z - gBeamSplitterZprime) / gBeamR;
689  }
690  }
691 }
static double gBeamR
static bool gFitBeamSplitters
static double gBeamSplitterZprime
static double gBSparam

◆ tecPlusFunction()

double TkLasBeamFitter::tecPlusFunction ( double *  x,
double *  par 
)
staticprivate

Definition at line 663 of file TkLasBeamFitter.cc.

References gBeamR, gBeamSplitterZprime, gBSparam, gFitBeamSplitters, x, and z.

Referenced by fitter().

663  {
664  double z = x[0]; // 'primed'? -> yes!!!
665 
666  if (z < gBeamSplitterZprime) {
667  return par[0] + par[1] * z;
668  } else {
669  if (gFitBeamSplitters) {
670  // par[2] = 2*tan(BeamSplitterAngle/2.0)
671  return par[0] + par[1] * z - par[2] * (z - gBeamSplitterZprime) / gBeamR;
672  } else {
673  return par[0] + par[1] * z - gBSparam * (z - gBeamSplitterZprime) / gBeamR;
674  }
675  }
676 }
static double gBeamR
static bool gFitBeamSplitters
static double gBeamSplitterZprime
static double gBSparam

◆ trackPhi()

void TkLasBeamFitter::trackPhi ( TkFittedLasBeam beam,
unsigned int &  hit,
double &  trackPhi,
double &  trackPhiRef,
double &  offset,
double &  slope,
double &  bsAngleParam,
double &  phiAtMinusParam,
double &  phiAtPlusParam,
double &  atThetaSplitParam,
std::vector< GlobalPoint > &  globHit 
)
private

Definition at line 812 of file TkLasBeamFitter.cc.

References funct::abs(), EcalCondDBWriter_cfi::beam, gBarrelModuleOffset, gBeamR, gBeamSplitterZprime, gBeamZ0, gHitsAtTecMinus, gHitZprime, gIsInnerBarrel, gTIBparam, gTOBparam, hltrates_dqm_sourceclient-live_cfg::offset, and slope.

Referenced by getLasBeams(), and globalTrackPoint().

822  {
823  // TECplus
824  if (beam.isTecInternal(1)) {
827  trackPhiRef = offset + slope * (gHitZprime[hit] + 1.0);
828  } else {
829  trackPhi = offset + slope * gHitZprime[hit] - bsAngleParam * (gHitZprime[hit] - gBeamSplitterZprime) / gBeamR;
830  trackPhiRef = offset + slope * (gHitZprime[hit] + 1.0) -
831  bsAngleParam * ((gHitZprime[hit] + 1.0) - gBeamSplitterZprime) / gBeamR;
832  }
833  }
834  // TECminus
835  else if (beam.isTecInternal(-1)) {
838  trackPhiRef = offset + slope * (gHitZprime[hit] + 1.0);
839  } else {
840  trackPhi = offset + slope * gHitZprime[hit] + bsAngleParam * (gHitZprime[hit] - gBeamSplitterZprime) / gBeamR;
841  trackPhiRef = offset + slope * (gHitZprime[hit] + 1.0) +
842  bsAngleParam * ((gHitZprime[hit] + 1.0) - gBeamSplitterZprime) / gBeamR;
843  }
844  }
845  // ATs
846  else {
847  // TECminus
848  if (gHitZprime[hit] < -gBeamSplitterZprime - 2.0 * gBeamZ0) {
850  trackPhiRef = offset + slope * (gHitZprime[hit] + 1.0);
851  }
852  // BarrelMinus
853  else if (-gBeamSplitterZprime - 2.0 * gBeamZ0 < gHitZprime[hit] && gHitZprime[hit] < -gBeamZ0) {
854  if (!gIsInnerBarrel) {
855  trackPhi = offset + slope * gHitZprime[hit] + gTOBparam * (phiAtMinusParam + atThetaSplitParam);
856  } else {
857  trackPhi = offset + slope * gHitZprime[hit] - gTIBparam * (phiAtMinusParam - atThetaSplitParam);
858  }
860  }
861  // BarrelPlus
863  if (!gIsInnerBarrel) {
864  trackPhi = offset + slope * gHitZprime[hit] + gTOBparam * (phiAtPlusParam - atThetaSplitParam);
865  } else {
866  trackPhi = offset + slope * gHitZprime[hit] - gTIBparam * (phiAtPlusParam + atThetaSplitParam);
867  }
869  }
870  // TECplus
871  else {
872  trackPhi = offset + slope * gHitZprime[hit] - bsAngleParam * (gHitZprime[hit] - gBeamSplitterZprime) / gBeamR;
873  trackPhiRef = offset + slope * (gHitZprime[hit] + 1.0) -
874  bsAngleParam * ((gHitZprime[hit] + 1.0) - gBeamSplitterZprime) / gBeamR;
875  }
876  }
877 }
static vector< double > gHitZprime
static float gTOBparam
static double gBeamR
static const double slope[3]
static bool gIsInnerBarrel
static double gBeamSplitterZprime
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
static unsigned int gHitsAtTecMinus
static vector< double > gBarrelModuleOffset
void trackPhi(TkFittedLasBeam &beam, unsigned int &hit, double &trackPhi, double &trackPhiRef, double &offset, double &slope, double &bsAngleParam, double &phiAtMinusParam, double &phiAtPlusParam, double &atThetaSplitParam, std::vector< GlobalPoint > &globHit)
static float gTIBparam
static double gBeamZ0

Member Data Documentation

◆ fieldHandle

edm::ESHandle<MagneticField> TkLasBeamFitter::fieldHandle
private

Definition at line 155 of file TkLasBeamFitter.cc.

Referenced by buildTrajectory(), and endRunProduce().

◆ fitBeamSplitters_

bool TkLasBeamFitter::fitBeamSplitters_
private

Definition at line 159 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), fitBeam(), fitter(), and getLasBeams().

◆ fs

edm::Service<TFileService> TkLasBeamFitter::fs
private

Definition at line 162 of file TkLasBeamFitter.cc.

Referenced by endRunProduce().

◆ gBarrelModuleOffset

vector< double > TkLasBeamFitter::gBarrelModuleOffset
staticprivate

Definition at line 167 of file TkLasBeamFitter.cc.

Referenced by fitBeam(), getLasBeams(), and trackPhi().

◆ gBarrelModuleRadius

vector< double > TkLasBeamFitter::gBarrelModuleRadius
staticprivate

Definition at line 166 of file TkLasBeamFitter.cc.

Referenced by getLasBeams().

◆ gBeamR

double TkLasBeamFitter::gBeamR = 0.0
staticprivate

◆ gBeamSplitterZprime

double TkLasBeamFitter::gBeamSplitterZprime = 0.0
staticprivate

◆ gBeamZ0

double TkLasBeamFitter::gBeamZ0 = 0.0
staticprivate

Definition at line 171 of file TkLasBeamFitter.cc.

Referenced by atFunction(), buildTrajectory(), fitBeam(), getLasBeams(), and trackPhi().

◆ gBSparam

double TkLasBeamFitter::gBSparam = 0.0
staticprivate

◆ geometry

edm::ESHandle<TrackerGeometry> TkLasBeamFitter::geometry
private

Definition at line 156 of file TkLasBeamFitter.cc.

◆ geomToken_

const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> TkLasBeamFitter::geomToken_
private

Definition at line 151 of file TkLasBeamFitter.cc.

Referenced by endRunProduce().

◆ gFitBeamSplitters

bool TkLasBeamFitter::gFitBeamSplitters = false
staticprivate

Definition at line 175 of file TkLasBeamFitter.cc.

Referenced by atFunction(), endRunProduce(), tecMinusFunction(), and tecPlusFunction().

◆ gHitsAtTecMinus

unsigned int TkLasBeamFitter::gHitsAtTecMinus = 0
staticprivate

Definition at line 173 of file TkLasBeamFitter.cc.

Referenced by fitBeam(), getLasBeams(), getLasHits(), globalTrackPoint(), and trackPhi().

◆ gHitZprime

vector< double > TkLasBeamFitter::gHitZprime
staticprivate

◆ gIsInnerBarrel

bool TkLasBeamFitter::gIsInnerBarrel = false
staticprivate

Definition at line 176 of file TkLasBeamFitter.cc.

Referenced by atFunction(), getLasBeams(), and trackPhi().

◆ gTIBparam

float TkLasBeamFitter::gTIBparam = 0.097614
staticprivate

Definition at line 168 of file TkLasBeamFitter.cc.

Referenced by atFunction(), and trackPhi().

◆ gTOBparam

float TkLasBeamFitter::gTOBparam = 0.034949
staticprivate

Definition at line 169 of file TkLasBeamFitter.cc.

Referenced by atFunction(), and trackPhi().

◆ h_bsAngle

TH1F* TkLasBeamFitter::h_bsAngle
private

Definition at line 179 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_bsAngleVsBeam

TH2F* TkLasBeamFitter::h_bsAngleVsBeam
private

Definition at line 181 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_chi2

TH1F * TkLasBeamFitter::h_chi2
private

Definition at line 179 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_chi2ndof

TH1F * TkLasBeamFitter::h_chi2ndof
private

Definition at line 179 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_hitX

TH1F * TkLasBeamFitter::h_hitX
private

Definition at line 179 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_hitXAt

TH1F * TkLasBeamFitter::h_hitXAt
private

Definition at line 179 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_hitXTecMinus

TH1F * TkLasBeamFitter::h_hitXTecMinus
private

Definition at line 179 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_hitXTecPlus

TH1F * TkLasBeamFitter::h_hitXTecPlus
private

Definition at line 179 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_hitXvsZAt

TH2F * TkLasBeamFitter::h_hitXvsZAt
private

Definition at line 181 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_hitXvsZTecMinus

TH2F * TkLasBeamFitter::h_hitXvsZTecMinus
private

Definition at line 181 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_hitXvsZTecPlus

TH2F * TkLasBeamFitter::h_hitXvsZTecPlus
private

Definition at line 181 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_pull

TH1F * TkLasBeamFitter::h_pull
private

Definition at line 179 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_res

TH1F * TkLasBeamFitter::h_res
private

Definition at line 179 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_resAt

TH1F * TkLasBeamFitter::h_resAt
private

Definition at line 179 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_resTecMinus

TH1F * TkLasBeamFitter::h_resTecMinus
private

Definition at line 179 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_resTecPlus

TH1F * TkLasBeamFitter::h_resTecPlus
private

Definition at line 179 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_resVsHitAt

TH2F * TkLasBeamFitter::h_resVsHitAt
private

Definition at line 181 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_resVsHitTecMinus

TH2F * TkLasBeamFitter::h_resVsHitTecMinus
private

Definition at line 181 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_resVsHitTecPlus

TH2F * TkLasBeamFitter::h_resVsHitTecPlus
private

Definition at line 181 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_resVsZAt

TH2F * TkLasBeamFitter::h_resVsZAt
private

Definition at line 181 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_resVsZTecMinus

TH2F * TkLasBeamFitter::h_resVsZTecMinus
private

Definition at line 181 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ h_resVsZTecPlus

TH2F * TkLasBeamFitter::h_resVsZTecPlus
private

Definition at line 181 of file TkLasBeamFitter.cc.

Referenced by endRunProduce(), and getLasBeams().

◆ laserBeams

edm::Handle<TkLasBeamCollection> TkLasBeamFitter::laserBeams
private

Definition at line 154 of file TkLasBeamFitter.cc.

Referenced by endRunProduce().

◆ magFieldToken_

const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> TkLasBeamFitter::magFieldToken_
private

Definition at line 150 of file TkLasBeamFitter.cc.

Referenced by endRunProduce().

◆ nAtParameters_

unsigned int TkLasBeamFitter::nAtParameters_
private

Definition at line 160 of file TkLasBeamFitter.cc.

Referenced by getLasBeams().

◆ src_

const edm::InputTag TkLasBeamFitter::src_
private

Definition at line 158 of file TkLasBeamFitter.cc.