diff --git a/Tracking/TrkFitter/TrkiPatFitter/TrkiPatFitter/MaterialAllocator.h b/Tracking/TrkFitter/TrkiPatFitter/TrkiPatFitter/MaterialAllocator.h index ee94bd319f7cea59d4ee6f8cf88b3374925ce845..a5505690b4780265615aa53ef06b6036d78c293f 100755 --- a/Tracking/TrkFitter/TrkiPatFitter/TrkiPatFitter/MaterialAllocator.h +++ b/Tracking/TrkFitter/TrkiPatFitter/TrkiPatFitter/MaterialAllocator.h @@ -62,12 +62,12 @@ public: FitParameters& fitParameters, Garbage_t& garbage) const; - // allocate material - void allocateMaterial (std::vector<FitMeasurement*>& measurements, - ParticleHypothesis particleHypothesis, - const FitParameters& fitParameters, - const TrackParameters& startParameters, - Garbage_t& garbage) const; + // allocate material + void allocateMaterial(std::vector<FitMeasurement*>& measurements, + ParticleHypothesis particleHypothesis, + FitParameters& fitParameters, + const TrackParameters& startParameters, + Garbage_t& garbage) const; // initialize scattering (needs to know X0 integral) void initializeScattering (std::vector<FitMeasurement*>& measurements) const; @@ -81,13 +81,13 @@ public: void orderMeasurements(std::vector<FitMeasurement*>& measurements, Amg::Vector3D startDirection, Amg::Vector3D startPosition) const; - - // has material been reallocated? - bool reallocateMaterial (std::vector<FitMeasurement*>& measurements, - const FitParameters& fitParameters, - Garbage_t& garbage) const; -private: + // has material been reallocated? + bool reallocateMaterial(std::vector<FitMeasurement*>& measurements, + FitParameters& fitParameters, + Garbage_t& garbage) const; + + private: // add material delimiters to control aggregation void addSpectrometerDelimiters (std::vector<FitMeasurement*>& measurements) const; @@ -106,11 +106,11 @@ private: Garbage_t& garbage) const; // allocate material in inner detector - void indetMaterial (std::vector<FitMeasurement*>& measurements, - ParticleHypothesis particleHypothesis, - const TrackParameters& startParameters, - Garbage_t& garbage) const; - + void indetMaterial(std::vector<FitMeasurement*>& measurements, + ParticleHypothesis particleHypothesis, + const TrackParameters& startParameters, + Garbage_t& garbage) const; + // material aggregation std::pair<FitMeasurement*,FitMeasurement*> materialAggregation ( const std::vector<const TrackStateOnSurface*>& material, @@ -129,11 +129,11 @@ private: void printMeasurements(std::vector<FitMeasurement*>& measurements) const; // allocate material in spectrometer - void spectrometerMaterial (std::vector<FitMeasurement*>& measurements, - ParticleHypothesis particleHypothesis, - const FitParameters& fitParameters, - const TrackParameters& startParameters, - Garbage_t& garbage) const; + void spectrometerMaterial(std::vector<FitMeasurement*>& measurements, + ParticleHypothesis particleHypothesis, + FitParameters& fitParameters, + const TrackParameters& startParameters, + Garbage_t& garbage) const; // Makes sure m_spectrometerEntrance is created, once only, and thread-safe void createSpectrometerEntranceOnce() const; diff --git a/Tracking/TrkFitter/TrkiPatFitter/src/MaterialAllocator.cxx b/Tracking/TrkFitter/TrkiPatFitter/src/MaterialAllocator.cxx index ce3c21a6eed06743426943cd89bb29e6b01152a9..44b907fe99bb27a6ed47de342d0597ca703d6559 100755 --- a/Tracking/TrkFitter/TrkiPatFitter/src/MaterialAllocator.cxx +++ b/Tracking/TrkFitter/TrkiPatFitter/src/MaterialAllocator.cxx @@ -99,14 +99,14 @@ namespace Trk // retrieve the necessary Extrapolators (muon tracking geometry is very picky!) ATH_CHECK( m_extrapolator.retrieve() ); ATH_MSG_DEBUG("Retrieved tool " << m_extrapolator); - + ATH_CHECK( m_intersector.retrieve() ); ATH_MSG_DEBUG("Retrieved tool " << m_intersector); // retrieve services ATH_CHECK( m_trackingGeometrySvc.retrieve() ); ATH_MSG_DEBUG("Retrieved Svc " << m_trackingGeometrySvc); - + // need to create the IndetExit and MuonEntrance TrackingVolumes ATH_CHECK( m_trackingVolumesSvc.retrieve() ); ATH_MSG_DEBUG("Retrieved Svc " << m_trackingVolumesSvc); @@ -605,7 +605,7 @@ namespace Trk void MaterialAllocator::allocateMaterial(std::vector<FitMeasurement*>& measurements, ParticleHypothesis particleHypothesis, - const FitParameters& fitParameters, + FitParameters& fitParameters, const TrackParameters& startParameters, Garbage_t& garbage) const { // different strategies used for indet and muon spectrometer @@ -705,9 +705,9 @@ namespace Trk // missing TrackingGeometrySvc - no leading material will be added m_messageHelper->printWarning(0); return nullptr; - } + } createSpectrometerEntranceOnce(); - + } // check input parameters are really in the spectrometer @@ -860,7 +860,7 @@ namespace Trk bool MaterialAllocator::reallocateMaterial(std::vector<FitMeasurement*>& measurements, - const FitParameters& parameters, + FitParameters& parameters, Garbage_t& garbage) const { ATH_MSG_DEBUG(" reallocateSpectrometerMaterial "); @@ -2189,7 +2189,7 @@ namespace Trk void MaterialAllocator::spectrometerMaterial(std::vector<FitMeasurement*>& measurements, ParticleHypothesis particleHypothesis, - const FitParameters& fitParameters, + FitParameters& fitParameters, const TrackParameters& startParameters, Garbage_t& garbage) const { // return if no MS measurement @@ -2302,9 +2302,9 @@ namespace Trk // missing TrackingGeometrySvc - no spectrometer material added m_messageHelper->printWarning(2); return; - } + } createSpectrometerEntranceOnce(); - + } // entranceParameters are at the MS entrance surface (0 if perigee downstream) diff --git a/Tracking/TrkFitter/TrkiPatFitterUtils/TrkiPatFitterUtils/FitParameters.h b/Tracking/TrkFitter/TrkiPatFitterUtils/TrkiPatFitterUtils/FitParameters.h index 2ab47c78b88b99986b20ff29c1111e33f99e4d9b..70698a484e2cfdb489a18cb0aa71c376d58743c7 100755 --- a/Tracking/TrkFitter/TrkiPatFitterUtils/TrkiPatFitterUtils/FitParameters.h +++ b/Tracking/TrkFitter/TrkiPatFitterUtils/TrkiPatFitterUtils/FitParameters.h @@ -105,19 +105,20 @@ public: void setPhiInstability (void); double sinPhi (void) const; double sinTheta (void) const; - Perigee* startingPerigee (void) const; - const TrackParameters* trackParameters (MsgStream& log, - const FitMeasurement& measurement, - bool withCovariance=false) const; - void update (const Amg::VectorX& differences); - void update (Amg::Vector3D position, - Amg::Vector3D direction, - double qOverP, - const Amg::MatrixX& leadingCovariance); - const Amg::Vector3D& vertex (void) const; - double z0 (void) const; - -private: + Perigee* startingPerigee(void) const; + //The following can update parameters + const TrackParameters* trackParameters(MsgStream& log, + const FitMeasurement& measurement, + bool withCovariance = false); + void update(const Amg::VectorX& differences); + void update(Amg::Vector3D position, + Amg::Vector3D direction, + double qOverP, + const Amg::MatrixX& leadingCovariance); + const Amg::Vector3D& vertex(void) const; + double z0(void) const; + + private: // assignment: no semantics, no implementation FitParameters &operator= (const FitParameters&); @@ -126,9 +127,9 @@ private: std::vector<double> m_alignmentOffset; std::vector<double> m_alignmentOffsetConstraint; double m_cosPhi; - mutable double m_cosPhi1; + double m_cosPhi1; double m_cosTheta; - mutable double m_cosTheta1; + double m_cosTheta1; double m_cotTheta; double m_d0; Amg::VectorX* m_differences; @@ -150,7 +151,7 @@ private: bool m_phiInstability; Amg::Vector3D m_position; double m_qOverP; - mutable double m_qOverP1; + double m_qOverP1; std::vector<double> m_scattererPhi; std::vector<double> m_scattererTheta; double m_sinPhi; diff --git a/Tracking/TrkFitter/TrkiPatFitterUtils/TrkiPatFitterUtils/FitProcedure.h b/Tracking/TrkFitter/TrkiPatFitterUtils/TrkiPatFitterUtils/FitProcedure.h index bd9a0c9b1e87cd231f5edd1535321fafbb803c7c..521dab0007785bdfbc45c94367f4c7118e853c21 100755 --- a/Tracking/TrkFitter/TrkiPatFitterUtils/TrkiPatFitterUtils/FitProcedure.h +++ b/Tracking/TrkFitter/TrkiPatFitterUtils/TrkiPatFitterUtils/FitProcedure.h @@ -66,7 +66,7 @@ public: // retrieve result Track* constructTrack (const std::vector<FitMeasurement*>& measurements, - const FitParameters& parameters, + FitParameters& parameters, const TrackInfo& trackInfo, const DataVector<const TrackStateOnSurface>* leadingTSOS = 0); diff --git a/Tracking/TrkFitter/TrkiPatFitterUtils/TrkiPatFitterUtils/IMaterialAllocator.h b/Tracking/TrkFitter/TrkiPatFitterUtils/TrkiPatFitterUtils/IMaterialAllocator.h index 2949072a0db5ef04d98616d85484a17a0918f727..e4320ade450d2034ed268284287e9b7b80d3f9b7 100755 --- a/Tracking/TrkFitter/TrkiPatFitterUtils/TrkiPatFitterUtils/IMaterialAllocator.h +++ b/Tracking/TrkFitter/TrkiPatFitterUtils/TrkiPatFitterUtils/IMaterialAllocator.h @@ -51,13 +51,13 @@ public: ParticleHypothesis particleHypothesis, FitParameters& fitParameters, Garbage_t& garbage) const = 0; - + /**IMaterialAllocator interface: allocate material */ - virtual void allocateMaterial (std::vector<FitMeasurement*>& measurements, - ParticleHypothesis particleHypothesis, - const FitParameters& fitParameters, - const TrackParameters& startParameters, - Garbage_t& garbage) const = 0; + virtual void allocateMaterial(std::vector<FitMeasurement*>& measurements, + ParticleHypothesis particleHypothesis, + FitParameters& fitParameters, + const TrackParameters& startParameters, + Garbage_t& garbage) const = 0; /**IMaterialAllocator interface: initialize scattering (needs to know X0 integral) */ virtual void initializeScattering (std::vector<FitMeasurement*>& measurements) const = 0; @@ -72,12 +72,11 @@ public: virtual void orderMeasurements (std::vector<FitMeasurement*>& measurements, Amg::Vector3D startDirection, Amg::Vector3D startPosition) const = 0; - - /**IMaterialAllocator interface: has material been reallocated? */ - virtual bool reallocateMaterial (std::vector<FitMeasurement*>& measurements, - const FitParameters& fitParameters, - Garbage_t& garbage) const = 0; + /**IMaterialAllocator interface: has material been reallocated? */ + virtual bool reallocateMaterial(std::vector<FitMeasurement*>& measurements, + FitParameters& fitParameters, + Garbage_t& garbage) const = 0; }; } // end of namespace diff --git a/Tracking/TrkFitter/TrkiPatFitterUtils/src/FitParameters.cxx b/Tracking/TrkFitter/TrkiPatFitterUtils/src/FitParameters.cxx index 68bcf57ca450bcf6afe73a6aaf5f697818c09d6e..626a509d8c02ba46ba7f980f7d7f4138536a7c3a 100755 --- a/Tracking/TrkFitter/TrkiPatFitterUtils/src/FitParameters.cxx +++ b/Tracking/TrkFitter/TrkiPatFitterUtils/src/FitParameters.cxx @@ -21,7 +21,7 @@ #include "TrkiPatFitterUtils/ParameterType.h" namespace Trk{ - + FitParameters::FitParameters (const Perigee& perigee) : m_cosPhi1 (0.), m_cosTheta1 (0.), @@ -49,7 +49,7 @@ FitParameters::FitParameters (const Perigee& perigee) m_vertex (perigee.associatedSurface().center()), m_z0 (perigee.position().z()) { - Amg::Vector3D momentum = perigee.momentum(); + Amg::Vector3D momentum = perigee.momentum(); double ptInv0 = 1./momentum.perp(); m_cosPhi = ptInv0*momentum.x(); m_sinPhi = ptInv0*momentum.y(); @@ -164,7 +164,7 @@ FitParameters::addAlignment (bool constrained, double angle, double offset) { m_alignmentAngleConstraint[m_numberAlignments] = angle; m_alignmentOffsetConstraint[m_numberAlignments] = offset; - } + } ++m_numberAlignments; } @@ -270,7 +270,7 @@ FitParameters::parameterDifference (const Amg::VectorX& parameters) const difference(0,4) = (parameters(5) - m_qOverP)*Gaudi::Units::TeV; return difference; } - + void FitParameters::performCutStep (double cutStep) { @@ -344,7 +344,7 @@ FitParameters::print (MsgStream& log) const << std::setw(11) << std::setprecision(6) << std::atan2(m_sinPhi,m_cosPhi) << " phi" << std::setw(11) << std::setprecision(6) << std::acos(m_cosTheta) - << " theta" + << " theta" << std::setw(13) << std::setprecision(3) << m_sinTheta/(m_qOverP*Gaudi::Units::GeV) << " pT (GeV)"; } @@ -394,11 +394,11 @@ FitParameters::printVerbose (MsgStream& log) const << std::setw(10) << std::setprecision(4) << differences(0) << " (0) " << std::setw(10) << std::setprecision(4) << differences(1) - << " (1) " + << " (1) " << std::setw(10) << std::setprecision(5) << differences(2) - << " (2) " + << " (2) " << std::setw(10) << std::setprecision(5) << differences(3) - << " (3) " + << " (3) " << std::setw(13) << std::setprecision(9) << differences(4)*Gaudi::Units::GeV/Gaudi::Units::TeV << " (4) "; if (m_fitEnergyDeposit) @@ -406,7 +406,7 @@ FitParameters::printVerbose (MsgStream& log) const << std::setw(13) << std::setprecision(9) << differences(5)*Gaudi::Units::GeV/Gaudi::Units::TeV << " (5) "; log << std::endl; - + if (m_numberAlignments) { log << " dAlign ==== "; @@ -453,17 +453,17 @@ FitParameters::printVerbose (MsgStream& log) const << std::setw(12) << std::setprecision(4) << m_d0 << " transverse impact " << std::setw(10) << std::setprecision(4) << m_position.z() - << " z0 " + << " z0 " << std::setw(10) << std::setprecision(6) << std::atan2(m_sinPhi,m_cosPhi) << std::setw(10) << std::setprecision(6) << m_cotTheta - << " phi,cotTheta " + << " phi,cotTheta " << std::setw(13) << std::setprecision(9) << m_qOverP/m_sinTheta << " inverse pT " << std::setw(12) << std::setprecision(6) << m_sinTheta/(m_qOverP*Gaudi::Units::GeV) << " signed pT "; if (m_fitEnergyDeposit) { - // TODO: should give fitted energy loss + // TODO: should give fitted energy loss log << std::endl << " E before/after energy deposit" << std::setw(12) << std::setprecision(3) << 1./std::abs(m_qOverP*Gaudi::Units::GeV) @@ -519,7 +519,7 @@ FitParameters::reset (const FitParameters& parameters) { // method is needed to complement copy in places where design uses // a reference to a FitParameter pointer - // essentially a copy, with history of previous iteration removed + // essentially a copy, with history of previous iteration removed m_cosPhi = parameters.m_cosPhi; m_cosPhi1 = parameters.m_cosPhi1; m_cosTheta = parameters.m_cosTheta; @@ -562,8 +562,8 @@ FitParameters::reset (const FitParameters& parameters) m_scattererPhi[s] = parameters.m_scattererPhi[s]; m_scattererTheta[s] = parameters.m_scattererTheta[s]; } - - // restore difference history + + // restore difference history delete m_differences; if (parameters.m_differences) { @@ -601,7 +601,7 @@ FitParameters::scatteringAngles (const FitMeasurement& fitMeasurement, int scatt scattererSigmaTheta); } } - + void FitParameters::setPhiInstability (void) { m_phiInstability = true; } @@ -614,7 +614,7 @@ FitParameters::startingPerigee (void) const double charge = 1.; if (m_qOverP < 0.) charge = -1.; Amg::Vector3D momentum(pT*m_cosPhi,pT*m_sinPhi,pT*m_cotTheta); - + return new Perigee(m_position, momentum, charge, @@ -624,7 +624,7 @@ FitParameters::startingPerigee (void) const const TrackParameters* FitParameters::trackParameters (MsgStream& log, const FitMeasurement& measurement, - bool withCovariance) const + bool withCovariance) { // make checks necessary for the TrackParameters to be computed // 1) a Surface is required @@ -640,7 +640,7 @@ FitParameters::trackParameters (MsgStream& log, log << MSG::WARNING << "FitParameters::trackParameters - invalid measurement" << endmsg; return nullptr; } - + // 3) the intersection position has to lie sufficiently close to the Surface const TrackSurfaceIntersection& intersection = measurement.intersection(FittedTrajectory); Amg::Vector2D localPos; @@ -651,7 +651,7 @@ FitParameters::trackParameters (MsgStream& log, log << MSG::WARNING << "FitParameters::trackParameters - globalToLocal failure" << endmsg; return nullptr; } - + // cache parameters at EnergyDeposit if (measurement.isEnergyDeposit()) { @@ -661,7 +661,7 @@ FitParameters::trackParameters (MsgStream& log, m_cosTheta1 = intersection.direction().z(); m_qOverP1 = measurement.qOverP(); } - + // propagate full covariance to form localCovariance AmgSymMatrix(5)* covMatrix = nullptr; if (withCovariance @@ -739,7 +739,7 @@ FitParameters::trackParameters (MsgStream& log, } theta = M_PI - theta; } - + // finally can create the appropriate 'concrete' TrackParameters const TrackParameters* parameters = nullptr; const StraightLineSurface* line = dynamic_cast<const StraightLineSurface*>(measurement.surface()); @@ -767,7 +767,7 @@ FitParameters::trackParameters (MsgStream& log, covMatrix); return parameters; } - + const CylinderSurface* cylinder = dynamic_cast<const CylinderSurface*>(measurement.surface()); if (cylinder) { @@ -793,7 +793,7 @@ FitParameters::trackParameters (MsgStream& log, covMatrix); return parameters; } - + const PerigeeSurface* peri = dynamic_cast<const PerigeeSurface*>(measurement.surface()); if (peri) { @@ -806,7 +806,7 @@ FitParameters::trackParameters (MsgStream& log, covMatrix); return parameters; } - + log << MSG::WARNING << "FitParameters::trackParameters - unrecognized surface" << endmsg; delete covMatrix; return nullptr; @@ -837,7 +837,7 @@ FitParameters::update (const Amg::VectorX& differences) (*a++) += differences(++align); (*o++) += differences(++align); } - + // scattering angles std::vector<double>::iterator p = m_scattererPhi.begin(); std::vector<double>::iterator t = m_scattererTheta.begin(); @@ -847,14 +847,14 @@ FitParameters::update (const Amg::VectorX& differences) (*p++) += differences(++scat); (*t++) += differences(++scat); } - + // qOverP, cotTheta if (m_fitMomentum) m_qOverP += differences(4)/Gaudi::Units::TeV; m_cotTheta -= differences(3)/(m_sinTheta*m_sinTheta); - + // impose charge conservation and decreasing energy if (m_fitEnergyDeposit) - { + { m_qOverP1 += differences(5)/Gaudi::Units::TeV; double deposit = 1./std::abs(m_qOverP) - 1./std::abs(m_qOverP1); if (std::abs(deposit) < std::abs(m_minEnergyDeposit) @@ -865,7 +865,7 @@ FitParameters::update (const Amg::VectorX& differences) if (m_qOverP1 < 0.) m_qOverP = -m_qOverP; } } - + // protect phi against some rounding instabilities double sinDPhi = differences(2); double cosDPhi = 0.; @@ -896,7 +896,7 @@ FitParameters::update (const Amg::VectorX& differences) m_vertex.y() + m_d0*m_cosPhi, m_z0); } - + void FitParameters::update (Amg::Vector3D position, Amg::Vector3D direction, diff --git a/Tracking/TrkFitter/TrkiPatFitterUtils/src/FitProcedure.cxx b/Tracking/TrkFitter/TrkiPatFitterUtils/src/FitProcedure.cxx index 7adacb872eb5e0fc2d9669d3c3a6955de46d106b..1c8f3332d7c3705f980948a0c7949ed2b9b2999d 100755 --- a/Tracking/TrkFitter/TrkiPatFitterUtils/src/FitProcedure.cxx +++ b/Tracking/TrkFitter/TrkiPatFitterUtils/src/FitProcedure.cxx @@ -41,7 +41,7 @@ #include "TrkiPatFitterUtils/MeasurementProcessor.h" namespace Trk{ - + // constructor FitProcedure::FitProcedure (bool constrainedAlignmentEffects, bool extendedDebug, @@ -109,7 +109,7 @@ FitProcedure::clear (void) Track* FitProcedure::constructTrack (const std::vector<FitMeasurement*>& measurements, - const FitParameters& parameters, + FitParameters& parameters, const TrackInfo& trackInfo, const DataVector<const TrackStateOnSurface>* leadingTSOS) { @@ -118,7 +118,7 @@ FitProcedure::constructTrack (const std::vector<FitMeasurement*>& measurements, // NB keep first and last measurements distinct i.e. separate TSOS (no scatterers etc) // NB trackParameters outwards from TSOS i.e. always last FitMeas on surface - + // create vector of TSOS - reserve upper limit for size (+1 as starts with perigee) DataVector<const TrackStateOnSurface>* trackStateOnSurfaces = new DataVector<const TrackStateOnSurface>; @@ -162,12 +162,12 @@ FitProcedure::constructTrack (const std::vector<FitMeasurement*>& measurements, } } } - + // then append the fitted TSOS for (auto *m : measurements) { if (m->isMaterialDelimiter()) continue; - + // push back previous TSOS when fresh surface reached if (m->surface() != surface || alignmentEffects || m->alignmentEffects()) { @@ -257,7 +257,7 @@ FitProcedure::constructTrack (const std::vector<FitMeasurement*>& measurements, typePattern = defaultPattern; alignmentEffects = nullptr; } - + measurementBase = m->measurementBase()->clone(); typePattern.set(TrackStateOnSurface::Measurement); if (m->isOutlier()) typePattern.set(TrackStateOnSurface::Outlier); @@ -285,11 +285,11 @@ FitProcedure::constructTrack (const std::vector<FitMeasurement*>& measurements, { const EnergyLoss* energyLoss = meot->energyLoss()->clone(); typeMaterial.set(Trk::MaterialEffectsBase::EnergyLossEffects); - if (m->numberDoF()) // fitted scatterer + if (m->numberDoF()) // fitted scatterer { materialEffects = new MaterialEffectsOnTrack(m->materialEffects()->thicknessInX0(), - parameters.scatteringAngles(*m,scatter), + parameters.scatteringAngles(*m,scatter), energyLoss, m->materialEffects()->associatedSurface(), typeMaterial); @@ -299,7 +299,7 @@ FitProcedure::constructTrack (const std::vector<FitMeasurement*>& measurements, { materialEffects = new MaterialEffectsOnTrack(m->materialEffects()->thicknessInX0(), - parameters.scatteringAngles(*m), + parameters.scatteringAngles(*m), energyLoss, m->materialEffects()->associatedSurface(), typeMaterial); @@ -311,7 +311,7 @@ FitProcedure::constructTrack (const std::vector<FitMeasurement*>& measurements, { materialEffects = new MaterialEffectsOnTrack(m->materialEffects()->thicknessInX0(), - parameters.scatteringAngles(*m,scatter), + parameters.scatteringAngles(*m,scatter), m->materialEffects()->associatedSurface(), typeMaterial); ++scatter; @@ -320,12 +320,12 @@ FitProcedure::constructTrack (const std::vector<FitMeasurement*>& measurements, { materialEffects = new MaterialEffectsOnTrack(m->materialEffects()->thicknessInX0(), - parameters.scatteringAngles(*m), + parameters.scatteringAngles(*m), m->materialEffects()->associatedSurface(), typeMaterial); } } - + typePattern.set(TrackStateOnSurface::Scatterer); } else @@ -341,14 +341,14 @@ FitProcedure::constructTrack (const std::vector<FitMeasurement*>& measurements, { typePattern.set(TrackStateOnSurface::Perigee); } - + // or alignment effects else if (m->alignmentEffects()) { const AlignmentEffectsOnTrack& AEOT = *m->alignmentEffects(); unsigned align = m->alignmentParameter() - 1; - *m_log << MSG::VERBOSE <<" Fitprocedure AEOT input deltaTranslation " << AEOT.deltaTranslation() << " deltaAngle " << AEOT.deltaAngle() << " output Trans " << parameters.alignmentOffset(align) << " deltaAngle " << parameters.alignmentAngle(align) << endmsg; + *m_log << MSG::VERBOSE <<" Fitprocedure AEOT input deltaTranslation " << AEOT.deltaTranslation() << " deltaAngle " << AEOT.deltaAngle() << " output Trans " << parameters.alignmentOffset(align) << " deltaAngle " << parameters.alignmentAngle(align) << endmsg; alignmentEffects = new Trk::AlignmentEffectsOnTrack(parameters.alignmentOffset(align), AEOT.sigmaDeltaTranslation(), @@ -358,7 +358,7 @@ FitProcedure::constructTrack (const std::vector<FitMeasurement*>& measurements, m->surface()); typePattern.set(TrackStateOnSurface::Alignment); } - + // passive types: hole for now else if (m->isPassive()) { @@ -419,7 +419,7 @@ FitProcedure::execute(bool asymmetricCaloEnergy, parameters->print(*m_log); *m_log << endmsg; } - + // choose appropriate intersector ToolHandle<IIntersector>& intersector = chooseIntersector(measurements,*parameters); @@ -457,13 +457,13 @@ FitProcedure::execute(bool asymmetricCaloEnergy, { m_fitMatrices->usePerigee(*measurements.front()); } - + // set requested options and initial values - double ptInvCut = 1./m_minPt; // protection against trapped particles + double ptInvCut = 1./m_minPt; // protection against trapped particles m_cutStep = true; m_convergence = false; m_nearConvergence = false; - + // keep best (original if not reasonable quality) results double bestChiSq = m_chiSqCut; FitParameters* bestParameters = nullptr; @@ -489,7 +489,7 @@ FitProcedure::execute(bool asymmetricCaloEnergy, if (m_extendedDebug) m_fitMatrices->checkPointers(*m_log); if (m_verbose) m_fitMatrices->printDerivativeMatrix(); } - + if (! m_fitMatrices->solveEquations()) { fitCode = 11; // fails matrix inversion @@ -514,8 +514,8 @@ FitProcedure::execute(bool asymmetricCaloEnergy, if (m_verbose && ! m_iteration) m_fitMatrices->printWeightMatrix(); } ++m_iteration; - - // report parameters + + // report parameters if (m_verbose) { *m_log << MSG::VERBOSE << " ===== start iteration " << m_iteration; @@ -529,7 +529,7 @@ FitProcedure::execute(bool asymmetricCaloEnergy, } parameters->printVerbose(*m_log); } - + // check for some error conditions (if none found yet) if (fitCode) { @@ -577,7 +577,7 @@ FitProcedure::execute(bool asymmetricCaloEnergy, m_fitMatrices->numberDoF(), parameters->numberScatterers(), m_worstMeasurement); - + if (m_debug) { if (m_verbose) *m_log << endmsg; @@ -588,7 +588,7 @@ FitProcedure::execute(bool asymmetricCaloEnergy, // have extrapolation and derivatives, calculate residual measurementProcessor.calculateResiduals(); - + // check for remaining error conditions. If OK then compute chisquared. if (m_iteration > m_maxIter && ! m_cutStep && for_iPatTrack) { @@ -596,7 +596,7 @@ FitProcedure::execute(bool asymmetricCaloEnergy, } else if (m_iteration == 4 && m_chiSq > 1000. && for_iPatTrack) { - fitCode = 7; // fail with too high chisquared + fitCode = 7; // fail with too high chisquared } else if (! fitCode) { @@ -629,7 +629,7 @@ FitProcedure::execute(bool asymmetricCaloEnergy, << m_iteration << ", numberOscillations " << parameters->numberOscillations() << endmsg; } - + // perform cutstep if (m_cutStep) { @@ -657,7 +657,7 @@ FitProcedure::execute(bool asymmetricCaloEnergy, bestParameters = new FitParameters(*parameters); parameters->resetOscillations(); } - + if (bestParameters && ((m_convergence && m_chiSq > bestChiSq + 0.5) || (parameters->phiInstability() && m_iteration == m_maxIter))) @@ -721,9 +721,9 @@ FitProcedure::execute(bool asymmetricCaloEnergy, m_chiSq = perigeeQuality->chiSquared() + m_chiSq * static_cast<double>(m_numberDoF); m_numberDoF += perigeeQuality->numberDoF(); - m_chiSq /= static_cast<double>(m_numberDoF); + m_chiSq /= static_cast<double>(m_numberDoF); } - + // probability of chisquared m_fitProbability = 1.; if (m_numberDoF > 0 && m_chiSq > 0.) @@ -768,7 +768,7 @@ FitProcedure::execute(bool asymmetricCaloEnergy, return *m_fitQuality; } - + Amg::MatrixX* FitProcedure::fullCovariance () const { @@ -776,7 +776,7 @@ FitProcedure::fullCovariance () const // return const_cast<Amg::MatrixX*>(m_fitMatrices->fullCovariance()); return nullptr; // NOT mig5 } - + void FitProcedure::setMinIterations (int minIter) { @@ -789,7 +789,7 @@ FitProcedure::calculateChiSq(std::vector<FitMeasurement*>& measurements) { // convergence criterion const double dChisqConv = 0.025; - + // compute total chisquared and sum of hit differences // flag hit with highest chisquared contribution (on entry if RoadFit) m_chiSq = 0.; @@ -803,10 +803,10 @@ FitProcedure::calculateChiSq(std::vector<FitMeasurement*>& measurements) // m_chiSq += m_fitMatrices->perigeeChiSquared(); // continue; // } - + double residual = m->residual(); double DiffSq = residual*residual; - m_chiSq += DiffSq; + m_chiSq += DiffSq; if (m->isPositionMeasurement()) { if (m->isDrift()) driftResidual += residual; @@ -843,7 +843,7 @@ FitProcedure::calculateChiSq(std::vector<FitMeasurement*>& measurements) m_chRatio2 = 0.; m_chiSqMin = m_chiSq; } - + m_chiSqOld = m_chiSqMin; double DChiSq = m_chiSqOld - m_chiSq; if (DChiSq > -dChisqConv) @@ -867,7 +867,7 @@ FitProcedure::calculateChiSq(std::vector<FitMeasurement*>& measurements) // fitlocal.dParMin[n] = fitlocal.dParam[n]; // } // } - + // if (m_cutTaken) // { // *m_log << " Cut ??? Chi2,DChiSq " << m_chiSq @@ -878,9 +878,9 @@ FitProcedure::calculateChiSq(std::vector<FitMeasurement*>& measurements) // << std::setiosflags(std::ios::fixed) // << " transverse impact parameter" // << std::setw(10) << std::setprecision(4) << parameters->d0() -// << " Z0" +// << " Z0" // << std::setw(10) << std::setprecision(4) << parameters->z0() -// << " 1/Pt" +// << " 1/Pt" // << std::setw(12) << std::setprecision(6) << parameters->ptInv0(); // // for (int n = 0; n < m_numberParameters; ++n) // // *m_log << " " << fitlocal.dParam[n]; @@ -889,19 +889,19 @@ FitProcedure::calculateChiSq(std::vector<FitMeasurement*>& measurements) // // *m_log << " " << fitlocal.dParMin[n]; // *m_log << std::endl; // } - + if (m_iteration > 0) { m_chRatio2 = m_chRatio1; - m_chRatio1 = m_chiSq/m_chiSqOld; + m_chRatio1 = m_chiSq/m_chiSqOld; } if (m_fitMatrices->numberDriftCircles()) { m_driftSumLast = m_driftSum; m_driftSum = driftResidual/static_cast<double>(m_fitMatrices->numberDriftCircles()); } - - // + + // // debugging info if (m_verbose) { @@ -915,7 +915,7 @@ FitProcedure::calculateChiSq(std::vector<FitMeasurement*>& measurements) << " ChSq Ratio1/2 " << std::setw(9) << std::setprecision(3) << m_chRatio1 << std::setw(10) << std::setprecision(3) << m_chRatio2 << std::endl << " driftResidual " << std::setw(9) << std::setprecision(3) << m_driftSum - << " #driftCircles " << m_fitMatrices->numberDriftCircles() << std::endl + << " #driftCircles " << m_fitMatrices->numberDriftCircles() << std::endl << " CutTaken " << m_cutTaken << std::endl << "----------------------------------" << std::endl << " "; @@ -940,9 +940,9 @@ FitProcedure::calculateChiSq(std::vector<FitMeasurement*>& measurements) } } } - - // - // check for possible convergence (nearConvergence forces extra iteration) + + // + // check for possible convergence (nearConvergence forces extra iteration) if (! m_cutStep && ! m_nCuts && ( m_chiSq < 0.1 @@ -955,7 +955,7 @@ FitProcedure::calculateChiSq(std::vector<FitMeasurement*>& measurements) m_convergence = true; } else - { + { m_nearConvergence = true; if (m_verbose) *m_log << MSG::VERBOSE << " near convergence " << endmsg; } @@ -964,7 +964,7 @@ FitProcedure::calculateChiSq(std::vector<FitMeasurement*>& measurements) { m_nearConvergence = false; } - + // else take cutstep if divergent or oscillating m_cutStep = false; // else if (m_nCuts < 2 || m_nCuts > 4) @@ -987,10 +987,10 @@ FitProcedure::calculateChiSq(std::vector<FitMeasurement*>& measurements) // if (m_verbose) // { -// *m_log << "----------------------------------" -// << " Debugging Info in ChiSquare method" -// << " cutstep - nCuts, -ve chisquared change " << m_nCuts -// << "----------------------------------" << std::endl; +// *m_log << "----------------------------------" +// << " Debugging Info in ChiSquare method" +// << " cutstep - nCuts, -ve chisquared change " << m_nCuts +// << "----------------------------------" << std::endl; // } // m_convergence = false; // m_chiSq = m_chiSqMin; @@ -1002,7 +1002,7 @@ FitProcedure::chooseIntersector (std::vector<FitMeasurement*>& measurements, const FitParameters& parameters) const { if (m_lineFit) return m_straightLineIntersector; - + // decide which intersector to use for curved tracks (default RungeKutta) // ToolHandle<IIntersector>& intersector = m_rungeKuttaIntersector; @@ -1012,7 +1012,7 @@ FitProcedure::chooseIntersector (std::vector<FitMeasurement*>& measurements, ++m) { if (! (**m).isPositionMeasurement()) continue; - if (! m_solenoidalIntersector->isValid(parameters.position(),(**m).position())) break; + if (! m_solenoidalIntersector->isValid(parameters.position(),(**m).position())) break; return m_solenoidalIntersector; } @@ -1034,7 +1034,7 @@ FitProcedure::reportQuality(const std::vector<FitMeasurement*>& measurements, { case 1: *m_log << " missing Trk::Surface "; - break; + break; case 2: *m_log << " too many measurements for fit matrix size: " << measurements.size(); @@ -1072,7 +1072,7 @@ FitProcedure::reportQuality(const std::vector<FitMeasurement*>& measurements, *m_log << " ill-defined cotTheta " << parameters.cotTheta() << " with difference " << parameters.difference(3) << " at iteration# "<< m_fitQuality->iterations(); - break; + break; case 11: *m_log << " singular matrix fails inversion:" << " at iteration# "<< m_fitQuality->iterations();