CMS 3D CMS Logo

Namespaces | Classes | Typedefs | Enumerations | Functions | Variables
mkfit Namespace Reference

Namespaces

 Config
 
 Const
 
 mini_propagators
 
 StdSeq
 

Classes

struct  axis
 
struct  axis_base
 
struct  axis_pow2
 
struct  axis_pow2_base
 
struct  axis_pow2_u1
 
struct  BeamSpot
 
struct  binnor
 
class  CandCloner
 
class  CcAlloc
 
class  CcPool
 
class  CombCandidate
 
class  ConfigJson
 
class  ConfigJsonPatcher
 
struct  DataFile
 
struct  DataFileHeader
 
struct  DeadRegion
 
class  Event
 
class  EventOfCombCandidates
 
class  EventOfHits
 
struct  ExecutionContext
 
class  FindingFoos
 
struct  FitVal
 
class  Hit
 
struct  HitMatch
 
struct  HitMatchPair
 
struct  HitOnTrack
 
struct  HoTNode
 
struct  IdxChi2List
 
class  IterationConfig
 
class  IterationLayerConfig
 
struct  IterationMaskIfc
 
struct  IterationMaskIfcBase
 
struct  IterationMaskIfcCmssw
 
class  IterationParams
 
class  IterationSeedPartition
 
class  IterationsInfo
 
struct  LayerControl
 
class  LayerInfo
 
class  LayerNumberConverter
 
class  LayerOfHits
 
class  MatriplexErrParPackerSlurpIn
 
class  MatriplexPackerSlurpIn
 
class  MatriplexTrackPackerPlexify
 
struct  MCHitInfo
 
struct  MeasurementState
 
class  MkBase
 
class  MkBuilder
 
class  MkBuilderWrapper
 
class  MkFinder
 
class  MkFitter
 
class  MkJob
 
struct  ModuleInfo
 
class  Pool
 
class  PropagationConfig
 
class  PropagationFlags
 
class  radix_sort
 
class  rectvec
 
struct  ReducedTrack
 
class  Shell
 
struct  sortTracksByPhiStruct
 
class  SteeringParams
 
class  Track
 
class  TrackBase
 
class  TrackCand
 
class  TrackerInfo
 
class  TrackExtra
 
struct  TrackState
 
class  TTreeValidation
 
struct  UpdateIndices
 
class  Validation
 
struct  WSR_Result
 

Typedefs

using clean_duplicates_cf = void(TrackVec &, const IterationConfig &)
 
using clean_duplicates_func = std::function< clean_duplicates_cf >
 
using clean_seeds_cf = int(TrackVec &, const IterationConfig &, const BeamSpot &)
 
using clean_seeds_func = std::function< clean_seeds_cf >
 
typedef std::map< std::string, std::pair< cleanOpts, std::string > > cleanOptsMap
 
typedef std::vector< DeadRegionDeadVec
 
typedef std::vector< EventEventVec
 
using filter_candidates_cf = bool(const TrackCand &, const MkJob &)
 
using filter_candidates_func = std::function< filter_candidates_cf >
 
typedef std::map< int, FitValFitValLayMap
 
typedef std::vector< std::pair< int, float > > FltLayerPairVec
 
typedef std::vector< int > HitIdxVec
 
typedef std::map< int, std::vector< int > > HitLayerMap
 
using HitVec = std::vector< Hit >
 
typedef std::vector< HitOnTrackHoTVec
 
typedef std::pair< int, float > idchi2Pair
 
typedef std::vector< idchi2Pairidchi2PairVec
 
typedef std::map< int, std::map< int, std::vector< int > > > LayIdxIDVecMapMap
 
typedef std::map< std::string, std::pair< matchOpts, std::string > > matchOptsMap
 
using MatriplexHitPacker = MatriplexErrParPackerSlurpIn< Hit, float >
 
using MatriplexHoTPacker = MatriplexPackerSlurpIn< HitOnTrack >
 
using MatriplexTrackPacker = MatriplexErrParPackerSlurpIn< TrackBase, float >
 
typedef std::vector< MCHitInfoMCHitInfoVec
 
typedef Matriplex::Matriplex< float, 2, 2, NNMPlex22
 
typedef Matriplex::Matriplex< float, 2, HH, NNMPlex2H
 
typedef Matriplex::MatriplexSym< float, 2, NNMPlex2S
 
typedef Matriplex::Matriplex< float, 2, 1, NNMPlex2V
 
typedef Matriplex::Matriplex< float, 5, 5, NNMPlex55
 
typedef Matriplex::Matriplex< float, 5, 6, NNMPlex56
 
typedef Matriplex::Matriplex< float, 6, 5, NNMPlex65
 
typedef Matriplex::Matriplex< float, HH, 2, NNMPlexH2
 
typedef Matriplex::Matriplex< float, HH, HH, NNMPlexHH
 
using MPlexHitIdx = Matriplex::Matriplex< int, MPlexHitIdxMax, 1, NN >
 
typedef Matriplex::Matriplex< float, HH, LL, NNMPlexHL
 
typedef Matriplex::MatriplexSym< float, HH, NNMPlexHS
 
typedef Matriplex::Matriplex< float, HH, 1, NNMPlexHV
 
typedef Matriplex::Matriplex< float, LL, 2, NNMPlexL2
 
typedef Matriplex::Matriplex< float, LL, HH, NNMPlexLH
 
typedef Matriplex::Matriplex< float, LL, LL, NNMPlexLL
 
typedef Matriplex::MatriplexSym< float, LL, NNMPlexLS
 
typedef Matriplex::Matriplex< float, LL, 1, NNMPlexLV
 
typedef Matriplex::Matriplex< bool, 1, 1, NNMPlexQB
 
typedef Matriplex::Matriplex< float, 1, 1, NNMPlexQF
 
typedef Matriplex::Matriplex< short, 1, 1, NNMPlexQH
 
using MPlexQHoT = Matriplex::Matriplex< HitOnTrack, 1, 1, NN >
 
typedef Matriplex::Matriplex< int, 1, 1, NNMPlexQI
 
typedef Matriplex::Matriplex< unsigned short, 1, 1, NNMPlexQUH
 
typedef Matriplex::Matriplex< unsigned int, 1, 1, NNMPlexQUI
 
typedef std::array< int, 2 > PairIdx
 
typedef std::vector< PairIdxPairIdxVec
 
using partition_seeds_cf = void(const TrackerInfo &, const TrackVec &, const EventOfHits &, IterationSeedPartition &)
 
using partition_seeds_func = std::function< partition_seeds_cf >
 
typedef std::vector< ReducedTrackRedTrackVec
 
typedef std::map< std::string, std::pair< seedOpts, std::string > > seedOptsMap
 
typedef std::pair< int, int > SimTkIDInfo
 
typedef ROOT::Math::SMatrix< float, 2 > SMatrix22
 
typedef ROOT::Math::SMatrix< float, 2, 6 > SMatrix26
 
typedef ROOT::Math::SMatrix< float, 3 > SMatrix33
 
typedef ROOT::Math::SMatrix< float, 3, 6 > SMatrix36
 
typedef ROOT::Math::SMatrix< float, 6, 2 > SMatrix62
 
typedef ROOT::Math::SMatrix< float, 6, 3 > SMatrix63
 
typedef ROOT::Math::SMatrix< float, 6 > SMatrix66
 
typedef ROOT::Math::SMatrix< float, 2, 2, ROOT::Math::MatRepSym< float, 2 > > SMatrixSym22
 
typedef ROOT::Math::SMatrix< float, 3, 3, ROOT::Math::MatRepSym< float, 3 > > SMatrixSym33
 
typedef ROOT::Math::SMatrix< float, 6, 6, ROOT::Math::MatRepSym< float, 6 > > SMatrixSym66
 
typedef ROOT::Math::SVector< float, 2 > SVector2
 
typedef ROOT::Math::SVector< float, 3 > SVector3
 
typedef ROOT::Math::SVector< float, 6 > SVector6
 
typedef std::unordered_map< int, FitValLayMapTkIDtoFitValLayMapMap
 
typedef std::unordered_map< int, int > TkIDToTkIDMap
 
typedef std::unordered_map< int, std::vector< int > > TkIDToTkIDVecMap
 
typedef std::unordered_map< int, TSLayerPairVecTkIDToTSLayerPairVecMap
 
typedef std::unordered_map< int, TrackStateTkIDToTSMap
 
typedef std::vector< TSVecTkIDToTSVecVec
 
using track_score_cf = float(const int nfoundhits, const int ntailholes, const int noverlaphits, const int nmisshits, const float chi2, const float pt, const bool inFindCandidates)
 
using track_score_func = std::function< track_score_cf >
 
typedef std::vector< TrackExtraTrackExtraVec
 
using TrackVec = std::vector< Track >
 
typedef std::vector< TrackVecTrackVecVec
 
typedef std::array< int, 3 > TripletIdx
 
typedef std::vector< TripletIdxTripletIdxVec
 
typedef std::map< int, std::unordered_set< int > > TrkIDLaySetMap
 
typedef std::vector< std::pair< int, TrackState > > TSLayerPairVec
 
typedef std::vector< TrackStateTSVec
 

Enumerations

enum  cleanOpts { noCleaning, cleanSeedsN2, cleanSeedsPure, cleanSeedsBadLabel }
 
enum  KalmanFilterOperation { KFO_Calculate_Chi2 = 1, KFO_Update_Params = 2, KFO_Local_Cov = 4 }
 
enum  matchOpts { trkParamBased, hitBased, labelBased }
 
enum  PropagationFlagsEnum { PF_none = 0, PF_use_param_b_field = 0x1, PF_apply_material = 0x2, PF_copy_input_state_on_fail = 0x4 }
 
enum  seedOpts { simSeeds, cmsswSeeds, findSeeds }
 
enum  TkLayout { TkLayout::phase0 = 0, TkLayout::phase1 = 1, TkLayout::phase2 = 2 }
 
enum  WithinSensitiveRegion_e {
  WSR_Undef = -1, WSR_Inside = 0, WSR_Edge, WSR_Outside,
  WSR_Failed
}
 

Functions

void applyMaterialEffects (const MPlexQF &hitsRl, const MPlexQF &hitsXi, const MPlexQF &propSign, const MPlexHV &plNrm, MPlexLS &outErr, MPlexLV &outPar, const int N_proc)
 
int calculateCharge (const Hit &hit0, const Hit &hit1, const Hit &hit2)
 
int calculateCharge (const float hit0_x, const float hit0_y, const float hit1_x, const float hit1_y, const float hit2_x, const float hit2_y)
 
float cdist (float a)
 
void CFMap (const MPlexHH &A, const MPlexHV &B, MPlexHV &C)
 
template<typename Vector , typename Matrix >
float computeHelixChi2 (const Vector &simV, const Vector &recoV, const Matrix &recoM, const bool diagOnly=false)
 
void conformalFitMPlex (bool fitting, MPlexQI seedID, MPlexLS &outErr, MPlexLV &outPar, const MPlexHV &msPar0, const MPlexHV &msPar1, const MPlexHV &msPar2)
 
template<typename Traits , typename HitCollection >
edm::ProductID convertHits (const Traits &traits, const HitCollection &hits, mkfit::HitVec &mkFitHits, std::vector< TrackingRecHit const *> &clusterIndexToHit, std::vector< float > &clusterChargeVec, const TrackerTopology &ttopo, const TransientTrackingRecHitBuilder &ttrhBuilder, const MkFitGeometry &mkFitGeom)
 
template<typename T >
T cube (T x)
 
template<typename Matrix >
void diagonalOnly (Matrix &m)
 
double dtime ()
 
template<typename Matrix >
void dumpMatrix (Matrix m)
 
void execTrackerInfoCreatorPlugin (const std::string &base, TrackerInfo &ti, IterationsInfo &ii, bool verbose)
 
void findSeedsByRoadSearch (TripletIdxConVec &seed_idcs, std::vector< LayerOfHits > &evt_lay_hits, int lay1_size, Event *&ev)
 
void from_json (const nlohmann::json &nlohmann_json_j, mkfit::LayerControl &nlohmann_json_t)
 
void from_json (const nlohmann::ordered_json &nlohmann_json_j, mkfit::LayerControl &nlohmann_json_t)
 
void from_json (const nlohmann::json &nlohmann_json_j, mkfit::SteeringParams &nlohmann_json_t)
 
void from_json (const nlohmann::ordered_json &nlohmann_json_j, mkfit::SteeringParams &nlohmann_json_t)
 
void from_json (const nlohmann::json &nlohmann_json_j, mkfit::IterationLayerConfig &nlohmann_json_t)
 
void from_json (const nlohmann::ordered_json &nlohmann_json_j, mkfit::IterationLayerConfig &nlohmann_json_t)
 
void from_json (const nlohmann::json &nlohmann_json_j, mkfit::IterationParams &nlohmann_json_t)
 
void from_json (const nlohmann::ordered_json &nlohmann_json_j, mkfit::IterationParams &nlohmann_json_t)
 
void from_json (const nlohmann::json &nlohmann_json_j, mkfit::IterationConfig &nlohmann_json_t)
 
void from_json (const nlohmann::ordered_json &nlohmann_json_j, mkfit::IterationConfig &nlohmann_json_t)
 
void from_json (const nlohmann::json &nlohmann_json_j, mkfit::IterationsInfo &nlohmann_json_t)
 
void from_json (const nlohmann::ordered_json &nlohmann_json_j, mkfit::IterationsInfo &nlohmann_json_t)
 
float getEta (float r, float z)
 
float getEta (float theta)
 
float getEta (float x, float y, float z)
 
float getEtaErr2 (float x, float y, float z, float exx, float eyy, float ezz, float exy, float exz, float eyz)
 
float getHypot (float x, float y)
 
float getInvRad2 (float x, float y)
 
float getInvRadErr2 (float x, float y, float exx, float eyy, float exy)
 
int getMatchBin (const float pt)
 
float getPhi (float x, float y)
 
float getPhiErr2 (float x, float y, float exx, float eyy, float exy)
 
float getPxPxErr2 (float ipt, float phi, float vipt, float vphi)
 
float getPyPyErr2 (float ipt, float phi, float vipt, float vphi)
 
float getPzPzErr2 (float ipt, float theta, float vipt, float vtheta)
 
float getRad2 (float x, float y)
 
float getRadErr2 (float x, float y, float exx, float eyy, float exy)
 
float getScoreCand (const track_score_func &score_func, const TrackCand &cand1, bool penalizeTailMissHits=false, bool inFindCandidates=false)
 
float getScoreCand (const track_score_func &score_func, const Track &cand1, bool penalizeTailMissHits=false, bool inFindCandidates=false)
 
float getScoreStruct (const track_score_func &score_func, const IdxChi2List &cand1)
 
float getScoreWorstPossible ()
 
float getTheta (float r, float z)
 
float getThetaErr2 (float x, float y, float z, float exx, float eyy, float ezz, float exy, float exz, float eyz)
 
void helixAtPlane (const MPlexLV &inPar, const MPlexQI &inChg, const MPlexHV &plPnt, const MPlexHV &plNrm, MPlexQF &pathL, MPlexLV &outPar, MPlexLL &errorProp, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &pflags)
 
void helixAtRFromIterativeCCS (const MPlexLV &inPar, const MPlexQI &inChg, const MPlexQF &msRad, MPlexLV &outPar, MPlexLL &errorProp, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &pflags)
 
void helixAtRFromIterativeCCSFullJac (const MPlexLV &inPar, const MPlexQI &inChg, const MPlexQF &msRad, MPlexLV &outPar, MPlexLL &errorProp, MPlexQI &outFailFlag, const int N_proc)
 
void helixAtRFromIterativeCCSFullJac (const MPlexLV &inPar, const MPlexQI &inChg, const MPlexQF &msRad, MPlexLV &outPar, MPlexLL &errorProp, const int N_proc)
 
void helixAtZ (const MPlexLV &inPar, const MPlexQI &inChg, const MPlexQF &msZ, MPlexLV &outPar, MPlexLL &errorProp, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &pflags)
 
float hipo (float x, float y)
 
float hipo_sqr (float x, float y)
 
void intersectThirdLayer (const float a, const float b, const float hit1_x, const float hit1_y, float &lay2_x, float &lay2_y)
 
constexpr bool isFinite (float x)
 
bool isStripQCompatible (int itrack, bool isBarrel, const MPlexLS &pErr, const MPlexLV &pPar, const MPlexHS &msErr, const MPlexHV &msPar)
 
void kalmanComputeChi2 (const MPlexLS &psErr, const MPlexLV &psPar, const MPlexQI &inChg, const MPlexHS &msErr, const MPlexHV &msPar, MPlexQF &outChi2, const int N_proc)
 
void kalmanComputeChi2Endcap (const MPlexLS &psErr, const MPlexLV &psPar, const MPlexQI &inChg, const MPlexHS &msErr, const MPlexHV &msPar, MPlexQF &outChi2, const int N_proc)
 
void kalmanComputeChi2Plane (const MPlexLS &psErr, const MPlexLV &psPar, const MPlexQI &inChg, const MPlexHS &msErr, const MPlexHV &msPar, const MPlexHV &plNrm, const MPlexHV &plDir, MPlexQF &outChi2, const int N_proc)
 
void kalmanOperation (const int kfOp, const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, MPlexLS &outErr, MPlexLV &outPar, MPlexQF &outChi2, const int N_proc)
 
void kalmanOperationEndcap (const int kfOp, const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, MPlexLS &outErr, MPlexLV &outPar, MPlexQF &outChi2, const int N_proc)
 
void kalmanOperationPlane (const int kfOp, const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, const MPlexHV &plNrm, const MPlexHV &plDir, MPlexLS &outErr, MPlexLV &outPar, MPlexQF &outChi2, const int N_proc)
 
void kalmanPropagateAndComputeChi2 (const MPlexLS &psErr, const MPlexLV &psPar, const MPlexQI &inChg, const MPlexHS &msErr, const MPlexHV &msPar, MPlexQF &outChi2, MPlexLV &propPar, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &propFlags, const bool propToHit)
 
void kalmanPropagateAndComputeChi2Endcap (const MPlexLS &psErr, const MPlexLV &psPar, const MPlexQI &inChg, const MPlexHS &msErr, const MPlexHV &msPar, MPlexQF &outChi2, MPlexLV &propPar, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &propFlags, const bool propToHit)
 
void kalmanPropagateAndComputeChi2Plane (const MPlexLS &psErr, const MPlexLV &psPar, const MPlexQI &inChg, const MPlexHS &msErr, const MPlexHV &msPar, const MPlexHV &plNrm, const MPlexHV &plDir, MPlexQF &outChi2, MPlexLV &propPar, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &propFlags, const bool propToHit)
 
void kalmanPropagateAndUpdate (const MPlexLS &psErr, const MPlexLV &psPar, MPlexQI &Chg, const MPlexHS &msErr, const MPlexHV &msPar, MPlexLS &outErr, MPlexLV &outPar, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &propFlags, const bool propToHit)
 
void kalmanPropagateAndUpdateEndcap (const MPlexLS &psErr, const MPlexLV &psPar, MPlexQI &Chg, const MPlexHS &msErr, const MPlexHV &msPar, MPlexLS &outErr, MPlexLV &outPar, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &propFlags, const bool propToHit)
 
void kalmanPropagateAndUpdatePlane (const MPlexLS &psErr, const MPlexLV &psPar, MPlexQI &Chg, const MPlexHS &msErr, const MPlexHV &msPar, const MPlexHV &plNrm, const MPlexHV &plDir, MPlexLS &outErr, MPlexLV &outPar, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &propFlags, const bool propToHit)
 
void kalmanUpdate (const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, MPlexLS &outErr, MPlexLV &outPar, const int N_proc)
 
void kalmanUpdateEndcap (const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, MPlexLS &outErr, MPlexLV &outPar, const int N_proc)
 
void kalmanUpdatePlane (const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, const MPlexHV &plNrm, const MPlexHV &plDir, MPlexLS &outErr, MPlexLV &outPar, const int N_proc)
 
template<class T , class U >
bool operator== (const CcAlloc< T > &a, const CcAlloc< U > &b)
 
bool passStripChargePCMfromTrack (int itrack, bool isBarrel, unsigned int pcm, unsigned int pcmMin, const MPlexLV &pPar, const MPlexHS &msErr)
 
void print (std::string_view label, const MeasurementState &s)
 
void print (const TrackState &s)
 
void print (std::string pfx, int itrack, const Track &trk, bool print_hits=false)
 
void print (std::string pfx, const TrackState &s)
 
void print (std::string pfx, int itrack, const Track &trk, const Event &ev)
 
void propagateHelixToPlaneMPlex (const MPlexLS &inErr, const MPlexLV &inPar, const MPlexQI &inChg, const MPlexHV &plPnt, const MPlexHV &plNrm, MPlexLS &outErr, MPlexLV &outPar, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &pflags, const MPlexQI *noMatEffPtr)
 
void propagateHelixToRMPlex (const MPlexLS &inErr, const MPlexLV &inPar, const MPlexQI &inChg, const MPlexQF &msRad, MPlexLS &outErr, MPlexLV &outPar, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &pflags, const MPlexQI *noMatEffPtr)
 
void propagateHelixToZMPlex (const MPlexLS &inErr, const MPlexLV &inPar, const MPlexQI &inChg, const MPlexQF &msZ, MPlexLS &outErr, MPlexLV &outPar, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &pflags, const MPlexQI *noMatEffPtr)
 
void propagateLineToRMPlex (const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, MPlexLS &outErr, MPlexLV &outPar, const int N_proc)
 
void run_OneIteration (const TrackerInfo &trackerInfo, const IterationConfig &itconf, const EventOfHits &eoh, const std::vector< const std::vector< bool > *> &hit_masks, MkBuilder &builder, TrackVec &seeds, TrackVec &out_tracks, bool do_seed_clean, bool do_backward_fit, bool do_remove_duplicates)
 
std::vector< double > runBtpCe_MultiIter (Event &ev, const EventOfHits &eoh, MkBuilder &builder, int n)
 
double runBuildingTestPlexBestHit (Event &ev, const EventOfHits &eoh, MkBuilder &builder)
 
double runBuildingTestPlexCloneEngine (Event &ev, const EventOfHits &eoh, MkBuilder &builder)
 
void runBuildingTestPlexDumbCMSSW (Event &ev, const EventOfHits &eoh, MkBuilder &builder)
 
double runBuildingTestPlexStandard (Event &ev, const EventOfHits &eoh, MkBuilder &builder)
 
double runFittingTestPlex (Event &ev, std::vector< Track > &rectracks)
 
void sincos4 (const float x, float &sin, float &cos)
 
bool sortByEta (const Hit &hit1, const Hit &hit2)
 
bool sortByHitsChi2 (const std::pair< Track, TrackState > &cand1, const std::pair< Track, TrackState > &cand2)
 
bool sortByHitsChi2 (const Track &cand1, const Track &cand2)
 
bool sortByPhi (const Hit &hit1, const Hit &hit2)
 
bool sortByScoreCand (const Track &cand1, const Track &cand2)
 
bool sortByScoreStruct (const IdxChi2List &cand1, const IdxChi2List &cand2)
 
bool sortByScoreTrackCand (const TrackCand &cand1, const TrackCand &cand2)
 
bool sortByZ (const Hit &hit1, const Hit &hit2)
 
bool sortIDsByChi2 (const idchi2Pair &cand1, const idchi2Pair &cand2)
 
bool sortTracksByEta (const Track &track1, const Track &track2)
 
bool sortTracksByPhi (const Track &track1, const Track &track2)
 
template<typename T >
T sqr (T x)
 
float squashPhiGeneral (float phi)
 
template<typename Vector >
void squashPhiGeneral (Vector &v)
 
float squashPhiMinimal (float phi)
 
void squashPhiMPlex (MPlexLV &par, const int N_proc)
 
void squashPhiMPlexGeneral (MPlexLV &par, const int N_proc)
 
void to_json (nlohmann::ordered_json &nlohmann_json_j, const mkfit::LayerControl &nlohmann_json_t)
 
void to_json (nlohmann::json &nlohmann_json_j, const mkfit::LayerControl &nlohmann_json_t)
 
void to_json (nlohmann::ordered_json &nlohmann_json_j, const mkfit::SteeringParams &nlohmann_json_t)
 
void to_json (nlohmann::json &nlohmann_json_j, const mkfit::SteeringParams &nlohmann_json_t)
 
void to_json (nlohmann::ordered_json &nlohmann_json_j, const mkfit::IterationLayerConfig &nlohmann_json_t)
 
void to_json (nlohmann::json &nlohmann_json_j, const mkfit::IterationLayerConfig &nlohmann_json_t)
 
void to_json (nlohmann::json &nlohmann_json_j, const mkfit::IterationParams &nlohmann_json_t)
 
void to_json (nlohmann::ordered_json &nlohmann_json_j, const mkfit::IterationParams &nlohmann_json_t)
 
void to_json (nlohmann::ordered_json &nlohmann_json_j, const mkfit::IterationConfig &nlohmann_json_t)
 
void to_json (nlohmann::json &nlohmann_json_j, const mkfit::IterationConfig &nlohmann_json_t)
 
void to_json (nlohmann::ordered_json &nlohmann_json_j, const mkfit::IterationsInfo &nlohmann_json_t)
 
void to_json (nlohmann::json &nlohmann_json_j, const mkfit::IterationsInfo &nlohmann_json_t)
 

Variables

bool g_debug = true
 
ExecutionContext g_exe_ctx
 
constexpr Matriplex::idx_t HH = 3
 
constexpr Matriplex::idx_t LL = 6
 
static constexpr int MPlexHitIdxMax = 16
 
constexpr Matriplex::idx_t NN = 8
 

Typedef Documentation

◆ clean_duplicates_cf

using mkfit::clean_duplicates_cf = typedef void(TrackVec &, const IterationConfig &)

Definition at line 30 of file FunctionTypes.h.

◆ clean_duplicates_func

using mkfit::clean_duplicates_func = typedef std::function<clean_duplicates_cf>

Definition at line 31 of file FunctionTypes.h.

◆ clean_seeds_cf

using mkfit::clean_seeds_cf = typedef int(TrackVec &, const IterationConfig &, const BeamSpot &)

Definition at line 21 of file FunctionTypes.h.

◆ clean_seeds_func

using mkfit::clean_seeds_func = typedef std::function<clean_seeds_cf>

Definition at line 22 of file FunctionTypes.h.

◆ cleanOptsMap

typedef std::map<std::string, std::pair<cleanOpts, std::string> > mkfit::cleanOptsMap

Definition at line 25 of file ConfigStandalone.h.

◆ DeadVec

typedef std::vector<DeadRegion> mkfit::DeadVec

Definition at line 280 of file Hit.h.

◆ EventVec

typedef std::vector<Event> mkfit::EventVec

Definition at line 84 of file Event.h.

◆ filter_candidates_cf

using mkfit::filter_candidates_cf = typedef bool(const TrackCand &, const MkJob &)

Definition at line 27 of file FunctionTypes.h.

◆ filter_candidates_func

using mkfit::filter_candidates_func = typedef std::function<filter_candidates_cf>

Definition at line 28 of file FunctionTypes.h.

◆ FitValLayMap

typedef std::map<int, FitVal> mkfit::FitValLayMap

Definition at line 24 of file TTreeValidation.h.

◆ FltLayerPairVec

typedef std::vector<std::pair<int, float> > mkfit::FltLayerPairVec

Definition at line 121 of file TrackExtra.h.

◆ HitIdxVec

typedef std::vector<int> mkfit::HitIdxVec

Definition at line 17 of file Track.h.

◆ HitLayerMap

typedef std::map<int, std::vector<int> > mkfit::HitLayerMap

Definition at line 18 of file Track.h.

◆ HitVec

typedef std::vector< Hit > mkfit::HitVec

Definition at line 11 of file MkFitHitWrapper.h.

◆ HoTVec

typedef std::vector<HitOnTrack> mkfit::HoTVec

Definition at line 272 of file Hit.h.

◆ idchi2Pair

typedef std::pair<int, float> mkfit::idchi2Pair

Definition at line 232 of file TrackExtra.cc.

◆ idchi2PairVec

typedef std::vector<idchi2Pair> mkfit::idchi2PairVec

Definition at line 233 of file TrackExtra.cc.

◆ LayIdxIDVecMapMap

typedef std::map<int, std::map<int, std::vector<int> > > mkfit::LayIdxIDVecMapMap

Definition at line 37 of file TrackExtra.h.

◆ matchOptsMap

typedef std::map<std::string, std::pair<matchOpts, std::string> > mkfit::matchOptsMap

Definition at line 29 of file ConfigStandalone.h.

◆ MatriplexHitPacker

Definition at line 136 of file MatriplexPackers.h.

◆ MatriplexHoTPacker

Definition at line 139 of file MatriplexPackers.h.

◆ MatriplexTrackPacker

Definition at line 137 of file MatriplexPackers.h.

◆ MCHitInfoVec

typedef std::vector<MCHitInfo> mkfit::MCHitInfoVec

Definition at line 119 of file Hit.h.

◆ MPlex22

typedef Matriplex::Matriplex<float, 2, 2, NN> mkfit::MPlex22

Definition at line 60 of file Matrix.h.

◆ MPlex2H

typedef Matriplex::Matriplex<float, 2, HH, NN> mkfit::MPlex2H

Definition at line 69 of file Matrix.h.

◆ MPlex2S

typedef Matriplex::MatriplexSym<float, 2, NN> mkfit::MPlex2S

Definition at line 62 of file Matrix.h.

◆ MPlex2V

typedef Matriplex::Matriplex<float, 2, 1, NN> mkfit::MPlex2V

Definition at line 61 of file Matrix.h.

◆ MPlex55

typedef Matriplex::Matriplex<float, 5, 5, NN> mkfit::MPlex55

Definition at line 56 of file Matrix.h.

◆ MPlex56

typedef Matriplex::Matriplex<float, 5, 6, NN> mkfit::MPlex56

Definition at line 57 of file Matrix.h.

◆ MPlex65

typedef Matriplex::Matriplex<float, 6, 5, NN> mkfit::MPlex65

Definition at line 58 of file Matrix.h.

◆ MPlexH2

typedef Matriplex::Matriplex<float, HH, 2, NN> mkfit::MPlexH2

Definition at line 68 of file Matrix.h.

◆ MPlexHH

typedef Matriplex::Matriplex<float, HH, HH, NN> mkfit::MPlexHH

Definition at line 52 of file Matrix.h.

◆ MPlexHitIdx

using mkfit::MPlexHitIdx = typedef Matriplex::Matriplex<int, MPlexHitIdxMax, 1, NN>

Definition at line 13 of file MkFitter.h.

◆ MPlexHL

typedef Matriplex::Matriplex<float, HH, LL, NN> mkfit::MPlexHL

Definition at line 65 of file Matrix.h.

◆ MPlexHS

typedef Matriplex::MatriplexSym<float, HH, NN> mkfit::MPlexHS

Definition at line 54 of file Matrix.h.

◆ MPlexHV

typedef Matriplex::Matriplex<float, HH, 1, NN> mkfit::MPlexHV

Definition at line 53 of file Matrix.h.

◆ MPlexL2

typedef Matriplex::Matriplex<float, LL, 2, NN> mkfit::MPlexL2

Definition at line 67 of file Matrix.h.

◆ MPlexLH

typedef Matriplex::Matriplex<float, LL, HH, NN> mkfit::MPlexLH

Definition at line 64 of file Matrix.h.

◆ MPlexLL

typedef Matriplex::Matriplex<float, LL, LL, NN> mkfit::MPlexLL

Definition at line 48 of file Matrix.h.

◆ MPlexLS

typedef Matriplex::MatriplexSym<float, LL, NN> mkfit::MPlexLS

Definition at line 50 of file Matrix.h.

◆ MPlexLV

typedef Matriplex::Matriplex<float, LL, 1, NN> mkfit::MPlexLV

Definition at line 49 of file Matrix.h.

◆ MPlexQB

typedef Matriplex::Matriplex<bool, 1, 1, NN> mkfit::MPlexQB

Definition at line 77 of file Matrix.h.

◆ MPlexQF

typedef Matriplex::Matriplex<float, 1, 1, NN> mkfit::MPlexQF

Definition at line 71 of file Matrix.h.

◆ MPlexQH

typedef Matriplex::Matriplex<short, 1, 1, NN> mkfit::MPlexQH

Definition at line 74 of file Matrix.h.

◆ MPlexQHoT

using mkfit::MPlexQHoT = typedef Matriplex::Matriplex<HitOnTrack, 1, 1, NN>

Definition at line 14 of file MkFitter.h.

◆ MPlexQI

typedef Matriplex::Matriplex<int, 1, 1, NN> mkfit::MPlexQI

Definition at line 72 of file Matrix.h.

◆ MPlexQUH

typedef Matriplex::Matriplex<unsigned short, 1, 1, NN> mkfit::MPlexQUH

Definition at line 75 of file Matrix.h.

◆ MPlexQUI

typedef Matriplex::Matriplex<unsigned int, 1, 1, NN> mkfit::MPlexQUI

Definition at line 73 of file Matrix.h.

◆ PairIdx

typedef std::array<int, 2> mkfit::PairIdx

Definition at line 39 of file TrackExtra.h.

◆ PairIdxVec

typedef std::vector<PairIdx> mkfit::PairIdxVec

Definition at line 40 of file TrackExtra.h.

◆ partition_seeds_cf

using mkfit::partition_seeds_cf = typedef void(const TrackerInfo &, const TrackVec &, const EventOfHits &, IterationSeedPartition &)

Definition at line 24 of file FunctionTypes.h.

◆ partition_seeds_func

using mkfit::partition_seeds_func = typedef std::function<partition_seeds_cf>

Definition at line 25 of file FunctionTypes.h.

◆ RedTrackVec

typedef std::vector<ReducedTrack> mkfit::RedTrackVec

Definition at line 35 of file TrackExtra.h.

◆ seedOptsMap

typedef std::map<std::string, std::pair<seedOpts, std::string> > mkfit::seedOptsMap

Definition at line 21 of file ConfigStandalone.h.

◆ SimTkIDInfo

typedef std::pair<int, int> mkfit::SimTkIDInfo

Definition at line 16 of file Track.h.

◆ SMatrix22

typedef ROOT::Math::SMatrix<float, 2> mkfit::SMatrix22

Definition at line 16 of file MatrixSTypes.h.

◆ SMatrix26

typedef ROOT::Math::SMatrix<float, 2, 6> mkfit::SMatrix26

Definition at line 23 of file MatrixSTypes.h.

◆ SMatrix33

typedef ROOT::Math::SMatrix<float, 3> mkfit::SMatrix33

Definition at line 12 of file MatrixSTypes.h.

◆ SMatrix36

typedef ROOT::Math::SMatrix<float, 3, 6> mkfit::SMatrix36

Definition at line 20 of file MatrixSTypes.h.

◆ SMatrix62

typedef ROOT::Math::SMatrix<float, 6, 2> mkfit::SMatrix62

Definition at line 24 of file MatrixSTypes.h.

◆ SMatrix63

typedef ROOT::Math::SMatrix<float, 6, 3> mkfit::SMatrix63

Definition at line 21 of file MatrixSTypes.h.

◆ SMatrix66

typedef ROOT::Math::SMatrix<float, 6> mkfit::SMatrix66

Definition at line 9 of file MatrixSTypes.h.

◆ SMatrixSym22

typedef ROOT::Math::SMatrix<float, 2, 2, ROOT::Math::MatRepSym<float, 2> > mkfit::SMatrixSym22

Definition at line 17 of file MatrixSTypes.h.

◆ SMatrixSym33

typedef ROOT::Math::SMatrix<float, 3, 3, ROOT::Math::MatRepSym<float, 3> > mkfit::SMatrixSym33

Definition at line 13 of file MatrixSTypes.h.

◆ SMatrixSym66

typedef ROOT::Math::SMatrix<float, 6, 6, ROOT::Math::MatRepSym<float, 6> > mkfit::SMatrixSym66

Definition at line 8 of file MatrixSTypes.h.

◆ SVector2

typedef ROOT::Math::SVector<float, 2> mkfit::SVector2

Definition at line 18 of file MatrixSTypes.h.

◆ SVector3

typedef ROOT::Math::SVector<float, 3> mkfit::SVector3

Definition at line 14 of file MatrixSTypes.h.

◆ SVector6

typedef ROOT::Math::SVector<float, 6> mkfit::SVector6

Definition at line 10 of file MatrixSTypes.h.

◆ TkIDtoFitValLayMapMap

typedef std::unordered_map<int, FitValLayMap> mkfit::TkIDtoFitValLayMapMap

Definition at line 25 of file TTreeValidation.h.

◆ TkIDToTkIDMap

typedef std::unordered_map<int, int> mkfit::TkIDToTkIDMap

Definition at line 124 of file TrackExtra.h.

◆ TkIDToTkIDVecMap

typedef std::unordered_map<int, std::vector<int> > mkfit::TkIDToTkIDVecMap

Definition at line 125 of file TrackExtra.h.

◆ TkIDToTSLayerPairVecMap

typedef std::unordered_map<int, TSLayerPairVec> mkfit::TkIDToTSLayerPairVecMap

Definition at line 127 of file TrackExtra.h.

◆ TkIDToTSMap

typedef std::unordered_map<int, TrackState> mkfit::TkIDToTSMap

Definition at line 126 of file TrackExtra.h.

◆ TkIDToTSVecVec

typedef std::vector<TSVec> mkfit::TkIDToTSVecVec

Definition at line 119 of file TrackExtra.h.

◆ track_score_cf

using mkfit::track_score_cf = typedef float(const int nfoundhits, const int ntailholes, const int noverlaphits, const int nmisshits, const float chi2, const float pt, const bool inFindCandidates)

Definition at line 39 of file FunctionTypes.h.

◆ track_score_func

using mkfit::track_score_func = typedef std::function<track_score_cf>

Definition at line 40 of file FunctionTypes.h.

◆ TrackExtraVec

typedef std::vector< TrackExtra > mkfit::TrackExtraVec

Definition at line 13 of file MkStandaloneSeqs.h.

◆ TrackVec

typedef std::vector< Track > mkfit::TrackVec

Definition at line 8 of file MkFitOutputWrapper.h.

◆ TrackVecVec

typedef std::vector<TrackVec> mkfit::TrackVecVec

Definition at line 597 of file Track.h.

◆ TripletIdx

typedef std::array<int, 3> mkfit::TripletIdx

Definition at line 41 of file TrackExtra.h.

◆ TripletIdxVec

typedef std::vector<TripletIdx> mkfit::TripletIdxVec

Definition at line 42 of file TrackExtra.h.

◆ TrkIDLaySetMap

typedef std::map<int, std::unordered_set<int> > mkfit::TrkIDLaySetMap

Definition at line 38 of file TrackExtra.h.

◆ TSLayerPairVec

typedef std::vector<std::pair<int, TrackState> > mkfit::TSLayerPairVec

Definition at line 120 of file TrackExtra.h.

◆ TSVec

typedef std::vector<TrackState> mkfit::TSVec

Definition at line 118 of file TrackExtra.h.

Enumeration Type Documentation

◆ cleanOpts

Enumerator
noCleaning 
cleanSeedsN2 
cleanSeedsPure 
cleanSeedsBadLabel 

Definition at line 24 of file ConfigStandalone.h.

◆ KalmanFilterOperation

Enumerator
KFO_Calculate_Chi2 
KFO_Update_Params 
KFO_Local_Cov 

Definition at line 11 of file KalmanUtilsMPlex.h.

◆ matchOpts

Enumerator
trkParamBased 
hitBased 
labelBased 

Definition at line 28 of file ConfigStandalone.h.

◆ PropagationFlagsEnum

Enumerator
PF_none 
PF_use_param_b_field 
PF_apply_material 
PF_copy_input_state_on_fail 

Definition at line 8 of file PropagationConfig.h.

◆ seedOpts

Enumerator
simSeeds 
cmsswSeeds 
findSeeds 

Definition at line 20 of file ConfigStandalone.h.

◆ TkLayout

enum mkfit::TkLayout
strong
Enumerator
phase0 
phase1 
phase2 

Definition at line 8 of file LayerNumberConverter.h.

◆ WithinSensitiveRegion_e

Enumerator
WSR_Undef 
WSR_Inside 
WSR_Edge 
WSR_Outside 
WSR_Failed 

Definition at line 17 of file TrackerInfo.h.

Function Documentation

◆ applyMaterialEffects()

void mkfit::applyMaterialEffects ( const MPlexQF hitsRl,
const MPlexQF hitsXi,
const MPlexQF propSign,
const MPlexHV plNrm,
MPlexLS outErr,
MPlexLV outPar,
const int  N_proc 
)

Definition at line 1239 of file PropagationMPlex.cc.

References funct::abs(), HLT_2023v12_cff::beta, funct::cos(), plot_hgcal_utils::dEdx, MillePedeFileConverter_cfg::e, f, CustomPhysics_cfi::gamma, Exhume::I, dqm-mbProfile::log, SiStripPI::max, hlt_dqm_clientPB-live_cfg::me, dqmiodumpmetadata::n, NN, AlCaHLTBitMon_ParallelJobs::p, SiStripOfflineCRack_cfg::p2, DiDispStaMuonMonitor_cfi::pt, funct::sin(), mathSSE::sqrt(), theta(), mkfit::Config::usePtMultScat, and cms::cuda::wmax.

Referenced by propagateHelixToPlaneMPlex(), propagateHelixToRMPlex(), and propagateHelixToZMPlex().

1245  {
1246 #pragma omp simd
1247  for (int n = 0; n < NN; ++n) {
1248  float radL = hitsRl.constAt(n, 0, 0);
1249  if (radL < 1e-13f)
1250  continue; //ugly, please fixme
1251  const float theta = outPar.constAt(n, 5, 0);
1252  // const float pt = 1.f / outPar.constAt(n, 3, 0); //fixme, make sure it is positive?
1253  const float ipt = outPar.constAt(n, 3, 0);
1254  const float pt = 1.f / ipt; //fixme, make sure it is positive?
1255  const float ipt2 = ipt * ipt;
1256  const float p = pt / std::sin(theta);
1257  const float pz = p * std::cos(theta);
1258  const float p2 = p * p;
1259  constexpr float mpi = 0.140; // m=140 MeV, pion
1260  constexpr float mpi2 = mpi * mpi; // m=140 MeV, pion
1261  const float beta2 = p2 / (p2 + mpi2);
1262  const float beta = std::sqrt(beta2);
1263  //radiation lenght, corrected for the crossing angle (cos alpha from dot product of radius vector and momentum)
1264  const float invCos =
1265  p / std::abs(pt * std::cos(outPar.constAt(n, 4, 0)) * plNrm.constAt(n, 0, 0) +
1266  pt * std::sin(outPar.constAt(n, 4, 0)) * plNrm.constAt(n, 1, 0) + pz * plNrm.constAt(n, 2, 0));
1267  radL = radL * invCos; //fixme works only for barrel geom
1268  // multiple scattering
1269  //vary independently phi and theta by the rms of the planar multiple scattering angle
1270  // XXX-KMD radL < 0, see your fixme above! Repeating bailout
1271  if (radL < 1e-13f)
1272  continue;
1273  // const float thetaMSC = 0.0136f*std::sqrt(radL)*(1.f+0.038f*std::log(radL))/(beta*p);// eq 32.15
1274  // const float thetaMSC2 = thetaMSC*thetaMSC;
1275  const float thetaMSC = 0.0136f * (1.f + 0.038f * std::log(radL)) / (beta * p); // eq 32.15
1276  const float thetaMSC2 = thetaMSC * thetaMSC * radL;
1277  if (Config::usePtMultScat) {
1278  outErr.At(n, 3, 3) += thetaMSC2 * pz * pz * ipt2 * ipt2;
1279  outErr.At(n, 3, 5) -= thetaMSC2 * pz * ipt2;
1280  outErr.At(n, 4, 4) += thetaMSC2 * p2 * ipt2;
1281  outErr.At(n, 5, 5) += thetaMSC2;
1282  } else {
1283  outErr.At(n, 4, 4) += thetaMSC2;
1284  outErr.At(n, 5, 5) += thetaMSC2;
1285  }
1286  //std::cout << "beta=" << beta << " p=" << p << std::endl;
1287  //std::cout << "multiple scattering thetaMSC=" << thetaMSC << " thetaMSC2=" << thetaMSC2 << " radL=" << radL << std::endl;
1288  // energy loss
1289  // XXX-KMD beta2 = 1 => 1 / sqrt(0)
1290  // const float gamma = 1.f/std::sqrt(1.f - std::min(beta2, 0.999999f));
1291  // const float gamma2 = gamma*gamma;
1292  const float gamma2 = (p2 + mpi2) / mpi2;
1293  const float gamma = std::sqrt(gamma2); //1.f/std::sqrt(1.f - std::min(beta2, 0.999999f));
1294  constexpr float me = 0.0005; // m=0.5 MeV, electron
1295  const float wmax = 2.f * me * beta2 * gamma2 / (1.f + 2.f * gamma * me / mpi + me * me / (mpi * mpi));
1296  constexpr float I = 16.0e-9 * 10.75;
1297  const float deltahalf = std::log(28.816e-9f * std::sqrt(2.33f * 0.498f) / I) + std::log(beta * gamma) - 0.5f;
1298  const float dEdx =
1299  beta < 1.f
1300  ? (2.f * (hitsXi.constAt(n, 0, 0) * invCos *
1301  (0.5f * std::log(2.f * me * beta2 * gamma2 * wmax / (I * I)) - beta2 - deltahalf) / beta2))
1302  : 0.f; //protect against infs and nans
1303  // dEdx = dEdx*2.;//xi in cmssw is defined with an extra factor 0.5 with respect to formula 27.1 in pdg
1304  //std::cout << "dEdx=" << dEdx << " delta=" << deltahalf << " wmax=" << wmax << " Xi=" << hitsXi.constAt(n,0,0) << std::endl;
1305  const float dP = propSign.constAt(n, 0, 0) * dEdx / beta;
1306  outPar.At(n, 3, 0) = p / (std::max(p + dP, 0.001f) * pt); //stay above 1MeV
1307  //assume 100% uncertainty
1308  outErr.At(n, 3, 3) += dP * dP / (p2 * pt * pt);
1309  }
1310  }
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
constexpr bool usePtMultScat
Definition: Config.h:52
T sqrt(T t)
Definition: SSEVec.h:19
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
const std::complex< double > I
Definition: I.h:8
double f[11][100]
Geom::Theta< T > theta() const
__host__ __device__ V V wmax

◆ calculateCharge() [1/2]

int mkfit::calculateCharge ( const Hit hit0,
const Hit hit1,
const Hit hit2 
)
inline

Definition at line 20 of file Track.h.

References mkfit::Hit::x(), and mkfit::Hit::y().

20  {
21  return ((hit2.y() - hit0.y()) * (hit2.x() - hit1.x()) > (hit2.y() - hit1.y()) * (hit2.x() - hit0.x()) ? 1 : -1);
22  }

◆ calculateCharge() [2/2]

int mkfit::calculateCharge ( const float  hit0_x,
const float  hit0_y,
const float  hit1_x,
const float  hit1_y,
const float  hit2_x,
const float  hit2_y 
)
inline

Definition at line 24 of file Track.h.

29  {
30  return ((hit2_y - hit0_y) * (hit2_x - hit1_x) > (hit2_y - hit1_y) * (hit2_x - hit0_x) ? 1 : -1);
31  }

◆ cdist()

float mkfit::cdist ( float  a)
inline

Definition at line 32 of file Config.h.

References a, mkfit::Const::PI, and mkfit::Const::TwoPI.

Referenced by mkfit::Event::clean_cms_seedtracks(), mkfit::StdSeq::clean_cms_seedtracks_iter(), mkfit::MkFinder::selectHitIndices(), and mkfit::MkFinder::selectHitIndicesV2().

32 { return a > Const::PI ? Const::TwoPI - a : a; }
constexpr float TwoPI
Definition: Config.h:8
constexpr float PI
Definition: Config.h:7
double a
Definition: hdecay.h:121

◆ CFMap()

void mkfit::CFMap ( const MPlexHH A,
const MPlexHV B,
MPlexHV C 
)
inline

Definition at line 89 of file ConformalUtilsMPlex.cc.

References A, a, ASSUME_ALIGNED, B, b, correctionTermsCaloMet_cff::C, HltBtagPostValidation_cff::c, N, and NN.

Referenced by conformalFitMPlex().

89  {
90  using idx_t = Matriplex::idx_t;
91 
92  // C = A * B, C is 3x1, A is 3x3 , B is 3x1
93 
94  typedef float T;
95  typedef float Tv;
96  const idx_t N = NN;
97 
98  const T* a = A.fArray;
99  ASSUME_ALIGNED(a, 64);
100  const Tv* b = B.fArray;
101  ASSUME_ALIGNED(b, 64);
102  Tv* c = C.fArray;
103  ASSUME_ALIGNED(c, 64);
104 
105 #include "RecoTracker/MkFitCore/standalone/CFMatrix33Vector3.ah"
106  }
Definition: APVGainStruct.h:7
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
#define N
Definition: blowfish.cc:9
double b
Definition: hdecay.h:120
double a
Definition: hdecay.h:121
Definition: APVGainStruct.h:7
long double T
#define ASSUME_ALIGNED(a, b)

◆ computeHelixChi2()

template<typename Vector , typename Matrix >
float mkfit::computeHelixChi2 ( const Vector simV,
const Vector recoV,
const Matrix &  recoM,
const bool  diagOnly = false 
)

Definition at line 652 of file Track.h.

References diagonalOnly(), and squashPhiGeneral().

Referenced by mkfit::TTreeValidation::fillEfficiencyTree(), mkfit::TTreeValidation::fillFakeRateTree(), mkfit::TrackExtra::setCMSSWTrackIDInfoByHits(), and mkfit::TrackExtra::setCMSSWTrackIDInfoByTrkParams().

652  {
653  Vector diffV = recoV - simV;
654  if (diffV.kSize > 2)
655  squashPhiGeneral(diffV);
656 
657  Matrix recoM_tmp = recoM;
658  if (diagOnly)
659  diagonalOnly(recoM_tmp);
660  int invFail(0);
661  const Matrix recoMI = recoM_tmp.InverseFast(invFail);
662 
663  return ROOT::Math::Dot(diffV * recoMI, diffV) / (diffV.kSize - 1);
664  }
ROOT::Math::Plane3D::Vector Vector
Definition: EcalHitMaker.cc:29
void squashPhiGeneral(Vector &v)
Definition: Track.h:645
portabletest::Matrix Matrix

◆ conformalFitMPlex()

void mkfit::conformalFitMPlex ( bool  fitting,
MPlexQI  seedID,
MPlexLS outErr,
MPlexLV outPar,
const MPlexHV msPar0,
const MPlexHV msPar1,
const MPlexHV msPar2 
)

Definition at line 111 of file ConformalUtilsMPlex.cc.

References A, a, funct::abs(), b, mkfit::Config::Bfield, correctionTermsCaloMet_cff::C, CFMap(), dcall, debug, dprint, dprintf, f, g_debug, getPhi(), getRad2(), getTheta(), hipo(), mps_fire::i, Matriplex::invertCramer(), dqmiolumiharvest::j, dqmdumpme::k, N, dqmiodumpmetadata::n, NN, mkfit::TrackState::parameters, mkfit::Config::phierr012, mkfit::Config::phierr049, mkfit::Const::PI3Over4, mkfit::Const::PIOver4, print(), pv::pT, mkfit::Config::ptinverr012, mkfit::Config::ptinverr049, multPhiCorr_741_25nsDY_cfi::px, multPhiCorr_741_25nsDY_cfi::py, diffTwoXMLs::r2, mkfit::Const::sol, mkfit::Config::thetaerr012, mkfit::Config::thetaerr049, findQualityFiles::v, mkfit::Config::varR, mkfit::Config::varXY, mkfit::Config::varZ, and x.

117  {
118  bool debug(false);
119 
120  using idx_t = Matriplex::idx_t;
121  const idx_t N = NN;
122 
123  // Store positions in mplex vectors... could consider storing in a 3x3 matrix, too
124  MPlexHV x, y, z, r2;
125 #pragma omp simd
126  for (int n = 0; n < N; ++n) {
127  x.At(n, 0, 0) = msPar0.constAt(n, 0, 0);
128  x.At(n, 1, 0) = msPar1.constAt(n, 0, 0);
129  x.At(n, 2, 0) = msPar2.constAt(n, 0, 0);
130 
131  y.At(n, 0, 0) = msPar0.constAt(n, 1, 0);
132  y.At(n, 1, 0) = msPar1.constAt(n, 1, 0);
133  y.At(n, 2, 0) = msPar2.constAt(n, 1, 0);
134 
135  z.At(n, 0, 0) = msPar0.constAt(n, 2, 0);
136  z.At(n, 1, 0) = msPar1.constAt(n, 2, 0);
137  z.At(n, 2, 0) = msPar2.constAt(n, 2, 0);
138 
139  for (int i = 0; i < 3; ++i) {
140  r2.At(n, i, 0) = getRad2(x.constAt(n, i, 0), y.constAt(n, i, 0));
141  }
142  }
143 
144  // Start setting the output parameters
145 #pragma omp simd
146  for (int n = 0; n < N; ++n) {
147  outPar.At(n, 0, 0) = x.constAt(n, 0, 0);
148  outPar.At(n, 1, 0) = y.constAt(n, 0, 0);
149  outPar.At(n, 2, 0) = z.constAt(n, 0, 0);
150  }
151 
152  // Use r-phi smearing to set initial error estimation for positions
153  // trackStates already initialized to identity for seeding ... don't store off-diag 0's, zero's for fitting set outside CF
154 #pragma omp simd
155  for (int n = 0; n < N; ++n) {
156  const float varPhi = Config::varXY / r2.constAt(n, 0, 0);
157  const float invvarR2 = Config::varR / r2.constAt(n, 0, 0);
158 
159  outErr.At(n, 0, 0) =
160  x.constAt(n, 0, 0) * x.constAt(n, 0, 0) * invvarR2 + y.constAt(n, 0, 0) * y.constAt(n, 0, 0) * varPhi;
161  outErr.At(n, 0, 1) = x.constAt(n, 0, 0) * y.constAt(n, 0, 0) * (invvarR2 - varPhi);
162 
163  outErr.At(n, 1, 0) = outErr.constAt(n, 0, 1);
164  outErr.At(n, 1, 1) =
165  y.constAt(n, 0, 0) * y.constAt(n, 0, 0) * invvarR2 + x.constAt(n, 0, 0) * x.constAt(n, 0, 0) * varPhi;
166 
167  outErr.At(n, 2, 2) = Config::varZ;
168  }
169 
170  MPlexQF initPhi;
171  MPlexQI xtou; // bool to determine "split space", i.e. map x to u or v
172 #pragma omp simd
173  for (int n = 0; n < N; ++n) {
174  initPhi.At(n, 0, 0) = std::abs(getPhi(x.constAt(n, 0, 0), y.constAt(n, 0, 0)));
175  xtou.At(n, 0, 0) =
176  ((initPhi.constAt(n, 0, 0) < Const::PIOver4 || initPhi.constAt(n, 0, 0) > Const::PI3Over4) ? 1 : 0);
177  }
178 
179  MPlexHV u, v;
180 #pragma omp simd
181  for (int n = 0; n < N; ++n) {
182  if (xtou.At(n, 0, 0)) // x mapped to u
183  {
184  for (int i = 0; i < 3; ++i) {
185  u.At(n, i, 0) = x.constAt(n, i, 0) / r2.constAt(n, i, 0);
186  v.At(n, i, 0) = y.constAt(n, i, 0) / r2.constAt(n, i, 0);
187  }
188  } else // x mapped to v
189  {
190  for (int i = 0; i < 3; ++i) {
191  u.At(n, i, 0) = y.constAt(n, i, 0) / r2.constAt(n, i, 0);
192  v.At(n, i, 0) = x.constAt(n, i, 0) / r2.constAt(n, i, 0);
193  }
194  }
195  }
196 
197  MPlexHH A;
198  //#pragma omp simd // triggers an internal compiler error with icc 18.0.2!
199  for (int n = 0; n < N; ++n) {
200  for (int i = 0; i < 3; ++i) {
201  A.At(n, i, 0) = 1.0f;
202  A.At(n, i, 1) = -u.constAt(n, i, 0);
203  A.At(n, i, 2) = -u.constAt(n, i, 0) * u.constAt(n, i, 0);
204  }
205  }
207  MPlexHV C;
208  CFMap(A, v, C);
209 
210  MPlexQF a, b;
211 #pragma omp simd
212  for (int n = 0; n < N; ++n) {
213  b.At(n, 0, 0) = 1.0f / (2.0f * C.constAt(n, 0, 0));
214  a.At(n, 0, 0) = b.constAt(n, 0, 0) * C.constAt(n, 1, 0);
215  }
216 
217  // constant used throughtout
218  const float k = (Const::sol * Config::Bfield) / 100.0f;
219 
220  MPlexQF vrx, vry, pT, px, py, pz;
221 #pragma omp simd
222  for (int n = 0; n < N; ++n) {
223  vrx.At(n, 0, 0) =
224  (xtou.constAt(n, 0, 0) ? x.constAt(n, 0, 0) - a.constAt(n, 0, 0) : x.constAt(n, 0, 0) - b.constAt(n, 0, 0));
225  vry.At(n, 0, 0) =
226  (xtou.constAt(n, 0, 0) ? y.constAt(n, 0, 0) - b.constAt(n, 0, 0) : y.constAt(n, 0, 0) - a.constAt(n, 0, 0));
227  pT.At(n, 0, 0) = k * hipo(vrx.constAt(n, 0, 0), vry.constAt(n, 0, 0));
228  px.At(n, 0, 0) = std::copysign(k * vry.constAt(n, 0, 0), x.constAt(n, 2, 0) - x.constAt(n, 0, 0));
229  py.At(n, 0, 0) = std::copysign(k * vrx.constAt(n, 0, 0), y.constAt(n, 2, 0) - y.constAt(n, 0, 0));
230  pz.At(n, 0, 0) = (pT.constAt(n, 0, 0) * (z.constAt(n, 2, 0) - z.constAt(n, 0, 0))) /
231  hipo((x.constAt(n, 2, 0) - x.constAt(n, 0, 0)), (y.constAt(n, 2, 0) - y.constAt(n, 0, 0)));
232  }
233 
234 #pragma omp simd
235  for (int n = 0; n < N; ++n) {
236  outPar.At(n, 3, 0) = 1.0f / pT.constAt(n, 0, 0);
237  outPar.At(n, 4, 0) = getPhi(px.constAt(n, 0, 0), py.constAt(n, 0, 0));
238  outPar.At(n, 5, 0) = getTheta(pT.constAt(n, 0, 0), pz.constAt(n, 0, 0));
239 #ifdef INWARDFIT // arctan is odd, so pz -> -pz means theta -> -theta
240  if (fitting)
241  outPar.At(n, 5, 0) *= -1.0f;
242 #endif
243  }
244 
245 #pragma omp simd
246  for (int n = 0; n < N; ++n) {
247  outErr.At(n, 3, 3) =
249  outErr.At(n, 4, 4) = (fitting ? Config::phierr049 * Config::phierr049 : Config::phierr012 * Config::phierr012);
250  outErr.At(n, 5, 5) =
252  }
253 
254  if (debug && g_debug) {
255  for (int n = 0; n < N; ++n) {
256  dprintf("afterCF seedID: %1u \n", seedID.constAt(n, 0, 0));
257  // do a dumb copy out
258  TrackState updatedState;
259  for (int i = 0; i < 6; i++) {
260  updatedState.parameters[i] = outPar.constAt(n, i, 0);
261  for (int j = 0; j < 6; j++) {
262  updatedState.errors[i][j] = outErr.constAt(n, i, j);
263  }
264  }
265 
266  dcall(print("CCS", updatedState));
267  updatedState.convertFromCCSToCartesian();
268  dcall(print("Pol", updatedState));
269  dprint("--------------------------------");
270  }
271  }
272  }
Matriplex::Matriplex< float, HH, HH, NN > MPlexHH
Definition: Matrix.h:52
Matriplex::Matriplex< float, HH, 1, NN > MPlexHV
Definition: Matrix.h:53
float getTheta(float r, float z)
Definition: Hit.h:36
constexpr float phierr012
float float float z
constexpr float thetaerr049
#define dcall(x)
Definition: Debug.h:97
bool g_debug
Definition: Debug.cc:2
constexpr float varZ
constexpr float varR
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
constexpr float Bfield
Definition: Config.h:55
float getRad2(float x, float y)
Definition: Hit.h:30
void print(TMatrixD &m, const char *label=nullptr, bool mathematicaFormat=false)
Definition: Utilities.cc:47
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
constexpr float ptinverr049
double f[11][100]
Matriplex::Matriplex< int, 1, 1, NN > MPlexQI
Definition: Matrix.h:72
constexpr float thetaerr012
#define debug
Definition: HDRShower.cc:19
constexpr float sol
Definition: Config.h:13
Matriplex::Matriplex< float, 1, 1, NN > MPlexQF
Definition: Matrix.h:71
#define N
Definition: blowfish.cc:9
void invertCramer(MPlex< T, D, D, N > &A, double *determ=nullptr)
Definition: Matriplex.h:730
float hipo(float x, float y)
Definition: Matrix.h:9
double b
Definition: hdecay.h:120
float getPhi(float x, float y)
Definition: Hit.h:34
#define dprint(x)
Definition: Debug.h:95
constexpr float varXY
double a
Definition: hdecay.h:121
constexpr float phierr049
void CFMap(const MPlexHH &A, const MPlexHV &B, MPlexHV &C)
float x
Definition: APVGainStruct.h:7
constexpr float PI3Over4
Definition: Config.h:11
constexpr float ptinverr012
#define dprintf(...)
Definition: Debug.h:98
constexpr float PIOver4
Definition: Config.h:10

◆ convertHits()

template<typename Traits , typename HitCollection >
edm::ProductID mkfit::convertHits ( const Traits &  traits,
const HitCollection hits,
mkfit::HitVec mkFitHits,
std::vector< TrackingRecHit const *> &  clusterIndexToHit,
std::vector< float > &  clusterChargeVec,
const TrackerTopology ttopo,
const TransientTrackingRecHitBuilder ttrhBuilder,
const MkFitGeometry mkFitGeom 
)

Definition at line 25 of file convertHits.h.

References ALCARECOTkAlJpsiMuMu_cff::charge, submitPVResolutionJobs::err, Exception, f, hfClusterShapes_cfi::hits, edm::ProductID::id(), TrackerTopology::isStereo(), TrackerTopology::layer(), LogTrace, MkFitGeometry::mkFitLayerNumber(), DetId::rawId(), DetId::subdetId(), MkFitGeometry::uniqueIdInLayer(), and UNLIKELY.

Referenced by MkFitSiPixelHitConverter::produce(), and MkFitSiStripHitConverter::produce().

32  {
33  if (hits.empty())
34  return edm::ProductID{};
35 
36  edm::ProductID clusterID;
37  {
38  const auto& lastClusterRef = hits.data().back().firstClusterRef();
39  clusterID = lastClusterRef.id();
40  if (lastClusterRef.index() >= mkFitHits.size()) {
41  auto const size = lastClusterRef.index();
42  mkFitHits.resize(size);
43  clusterIndexToHit.resize(size, nullptr);
44  if constexpr (Traits::applyCCC()) {
45  clusterChargeVec.resize(size, -1.f);
46  }
47  }
48  }
49 
50  for (const auto& detset : hits) {
51  const DetId detid = detset.detId();
52  const auto ilay = mkFitGeom.mkFitLayerNumber(detid);
53 
54  for (const auto& hit : detset) {
55  const auto charge = traits.clusterCharge(hit, detid);
56  if (!traits.passCCC(charge))
57  continue;
58 
59  const auto& gpos = hit.globalPosition();
60  SVector3 pos(gpos.x(), gpos.y(), gpos.z());
61  const auto& gerr = hit.globalPositionError();
63  err.At(0, 0) = gerr.cxx();
64  err.At(1, 1) = gerr.cyy();
65  err.At(2, 2) = gerr.czz();
66  err.At(0, 1) = gerr.cyx();
67  err.At(0, 2) = gerr.czx();
68  err.At(1, 2) = gerr.czy();
69 
70  auto clusterRef = hit.firstClusterRef();
71  if UNLIKELY (clusterRef.id() != clusterID) {
72  throw cms::Exception("LogicError")
73  << "Input hit collection has Refs to many cluster collections. Last hit had Ref to product " << clusterID
74  << ", but encountered Ref to product " << clusterRef.id() << " on detid " << detid.rawId();
75  }
76  const auto clusterIndex = clusterRef.index();
77  LogTrace("MkFitHitConverter") << "Adding hit detid " << detid.rawId() << " subdet " << detid.subdetId()
78  << " layer " << ttopo.layer(detid) << " isStereo " << ttopo.isStereo(detid)
79  << " zplus "
80  << " index " << clusterIndex << " ilay " << ilay;
81 
82  if UNLIKELY (clusterIndex >= mkFitHits.size()) {
83  mkFitHits.resize(clusterIndex + 1);
84  clusterIndexToHit.resize(clusterIndex + 1, nullptr);
85  if constexpr (Traits::applyCCC()) {
86  clusterChargeVec.resize(clusterIndex + 1, -1.f);
87  }
88  }
89  mkFitHits[clusterIndex] = mkfit::Hit(pos, err);
90  clusterIndexToHit[clusterIndex] = &hit;
91  if constexpr (Traits::applyCCC()) {
92  clusterChargeVec[clusterIndex] = charge;
93  }
94 
95  const auto uniqueIdInLayer = mkFitGeom.uniqueIdInLayer(ilay, detid.rawId());
96  traits.setDetails(mkFitHits[clusterIndex], *(hit.cluster()), uniqueIdInLayer, charge);
97  }
98  }
99  return clusterID;
100  }
size
Write out results.
unsigned int uniqueIdInLayer(int layer, unsigned int detId) const
Definition: MkFitGeometry.h:36
bool isStereo(const DetId &id) const
unsigned int layer(const DetId &id) const
#define LogTrace(id)
double f[11][100]
constexpr int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:48
SeedingHitSet::ConstRecHitPointer Hit
ProductIndex id() const
Definition: ProductID.h:35
Definition: DetId.h:17
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:57
int mkFitLayerNumber(DetId detId) const
#define UNLIKELY(x)
Definition: Likely.h:21
ROOT::Math::SVector< double, 3 > SVector3
Definition: V0Fitter.cc:47
ROOT::Math::SMatrix< float, 3, 3, ROOT::Math::MatRepSym< float, 3 > > SMatrixSym33
Definition: MatrixSTypes.h:13

◆ cube()

template<typename T >
T mkfit::cube ( T  x)
inline

Definition at line 18 of file Hit.h.

References x.

Referenced by getInvRadErr2().

18  {
19  return x * x * x;
20  }
float x

◆ diagonalOnly()

template<typename Matrix >
void mkfit::diagonalOnly ( Matrix &  m)
inline

Definition at line 27 of file MatrixSTypes.h.

References HltBtagPostValidation_cff::c, and visualization-live-secondInstance_cfg::m.

Referenced by computeHelixChi2().

27  {
28  for (int r = 0; r < m.kRows; r++) {
29  for (int c = 0; c < m.kCols; c++) {
30  if (r != c)
31  m[r][c] = 0.f;
32  }
33  }
34  }

◆ dtime()

double mkfit::dtime ( )
inline

◆ dumpMatrix()

template<typename Matrix >
void mkfit::dumpMatrix ( Matrix  m)

Definition at line 37 of file MatrixSTypes.h.

References HltBtagPostValidation_cff::c, gather_cfg::cout, and visualization-live-secondInstance_cfg::m.

Referenced by print().

37  {
38  for (int r = 0; r < m.kRows; ++r) {
39  for (int c = 0; c < m.kCols; ++c) {
40  std::cout << std::setw(12) << m.At(r, c) << " ";
41  }
42  std::cout << std::endl;
43  }
44  }

◆ execTrackerInfoCreatorPlugin()

void mkfit::execTrackerInfoCreatorPlugin ( const std::string &  base,
TrackerInfo ti,
IterationsInfo ii,
bool  verbose 
)

Definition at line 92 of file ConfigStandalone.cc.

References newFWLiteAna::base, beamvalidation::exit(), custom_jme_cff::foo, h, cuy::ii, castor_dqm_sourceclient_file_cfg::path, mkfit::TrackerInfo::read_bin_file(), edm_modernize_messagelogger::stat, submitPVResolutionJobs::stderr, AlCaHLTBitMon_QueryRunRegistry::string, and verbose.

Referenced by initGeom().

92  {
93  const std::string soname = base + ".so";
94  const std::string binname = base + ".bin";
95  struct stat st;
96 
97  int si = 0;
98  while (search_path[si]) {
100  const char *envpath = std::getenv("MKFIT_BASE");
101  if (envpath != nullptr) {
102  path += envpath;
103  path += "/";
104  }
105  path += search_path[si];
106  std::string sopath = path + soname;
107  if (stat(sopath.c_str(), &st) == 0) {
108  printf("execTrackerInfoCreatorPlugin processing '%s'\n", sopath.c_str());
109 
110  void *h = dlopen(sopath.c_str(), RTLD_LAZY);
111  if (!h) {
112  perror("dlopen failed");
113  exit(2);
114  }
115 
116  long long *p2f = (long long *)dlsym(h, "TrackerInfoCreator_ptr");
117  if (!p2f) {
118  perror("dlsym failed");
119  exit(2);
120  }
121 
122  std::string binpath = path + binname;
123  int binsr = stat(binpath.c_str(), &st);
124  printf("execTrackerInfoCreatorPlugin has%s found TrackerInfo binary file '%s'\n",
125  binsr ? " NOT" : "",
126  binpath.c_str());
127  if (binsr == 0)
128  ti.read_bin_file(binpath);
129 
130  TrackerInfoCreator_foo foo = (TrackerInfoCreator_foo)(*p2f);
131  foo(ti, ii, verbose);
132 
133  return;
134  }
135 
136  ++si;
137  }
138 
139  fprintf(stderr, "TrackerInfo plugin '%s' not found in search path.\n", soname.c_str());
140  exit(2);
141  }
bool verbose
base
Main Program
Definition: newFWLiteAna.py:92
ii
Definition: cuy.py:589
The Signals That Services Can Subscribe To This is based on ActivityRegistry h
Helper function to determine trigger accepts.
Definition: Activities.doc:4
def exit(msg="")

◆ findSeedsByRoadSearch()

void mkfit::findSeedsByRoadSearch ( TripletIdxConVec &  seed_idcs,
std::vector< LayerOfHits > &  evt_lay_hits,
int  lay1_size,
Event *&  ev 
)

Definition at line 34 of file seedtestMPlex.cc.

References a, funct::abs(), b, debug, mkfit::Config::dEtaSeedTrip, dprint, makeMEIFBenchmarkPlots::ev, JetChargeProducer_cfi::exp, f, getHypot(), getPhi(), getRad2(), mps_fire::i, heavyIonCSV_trainingSettings::idx, intersectThirdLayer(), mkfit::Config::lay01angdiff, SiStripPI::max, mkfit::Config::maxCurvR, mkfit::Hit::mcTrackID(), HLTObjectMonitor_cfi::mr, TtSemiLepEvtBuilder_cfi::mt, mkfit::Config::numHitsPerTask, mkfit::Hit::phi(), mkfit::LayerOfHits::refHit(), mkfit::Config::seed_d0cut, mkfit::Config::seed_z0cut, mkfit::Config::seed_z1cut, mathSSE::sqrt(), funct::tan(), mkfit::Hit::x(), mkfit::Hit::y(), and mkfit::Hit::z().

37  {
38 #ifdef DEBUG
39  bool debug(false);
40 #endif
41 
42  // MIMI hack: Config::nlayers_per_seed = 4
43  // const float seed_z2cut = (Config::nlayers_per_seed * Config::fRadialSpacing) / std::tan(2.0f*std::atan(std::exp(-1.0f*Config::dEtaSeedTrip)));
44 #ifdef DEBUG
45  const float seed_z2cut =
46  (4 * Config::fRadialSpacing) / std::tan(2.0f * std::atan(std::exp(-1.0f * Config::dEtaSeedTrip)));
47 #endif
48 
49  // 0 = first layer, 1 = second layer, 2 = third layer
50  const LayerOfHits& lay1_hits = evt_lay_hits[1];
51  LayerOfHits& lay0_hits = evt_lay_hits[0];
52  LayerOfHits& lay2_hits = evt_lay_hits[2];
53 
54  tbb::parallel_for(
55  tbb::blocked_range<int>(0, lay1_size, std::max(1, Config::numHitsPerTask)),
56  [&](const tbb::blocked_range<int>& i) {
57  TripletIdxVec temp_thr_seed_idcs;
58  for (int ihit1 = i.begin(); ihit1 < i.end(); ++ihit1) {
59  const Hit& hit1 = lay1_hits.refHit(ihit1);
60  const float hit1_z = hit1.z();
61 
62  dprint("ihit1: " << ihit1 << " mcTrackID: " << hit1.mcTrackID(ev->simHitsInfo_) << " phi: " << hit1.phi()
63  << " z: " << hit1.z());
64  dprint(" predphi: " << hit1.phi() << "+/-" << Config::lay01angdiff << " predz: " << hit1.z() / 2.0f << "+/-"
65  << Config::seed_z0cut / 2.0f << std::endl);
66 
67  std::vector<int> cand_hit0_indices; // pass by reference
68  // MIMI lay0_hits.selectHitIndices(hit1_z/2.0f,hit1.phi(),Config::seed_z0cut/2.0f,Config::lay01angdiff,cand_hit0_indices,true,false);
69  // loop over first layer hits
70  for (auto&& ihit0 : cand_hit0_indices) {
71  const Hit& hit0 = lay0_hits.refHit(ihit0);
72  const float hit0_z = hit0.z();
73  const float hit0_x = hit0.x();
74  const float hit0_y = hit0.y();
75  const float hit1_x = hit1.x();
76  const float hit1_y = hit1.y();
77  const float hit01_r2 = getRad2(hit0_x - hit1_x, hit0_y - hit1_y);
78 
79  const float quad = std::sqrt((4.0f * Config::maxCurvR * Config::maxCurvR - hit01_r2) / hit01_r2);
80 
81  // center of negative curved track
82  const float aneg = 0.5f * ((hit0_x + hit1_x) - (hit0_y - hit1_y) * quad);
83  const float bneg = 0.5f * ((hit0_y + hit1_y) + (hit0_x - hit1_x) * quad);
84 
85  // negative points of intersection with third layer
86  float lay2_negx = 0.0f, lay2_negy = 0.0f;
87  intersectThirdLayer(aneg, bneg, hit1_x, hit1_y, lay2_negx, lay2_negy);
88 #ifdef DEBUG
89  const float lay2_negphi = getPhi(lay2_negx, lay2_negy);
90 #endif
91 
92  // center of positive curved track
93  const float apos = 0.5f * ((hit0_x + hit1_x) + (hit0_y - hit1_y) * quad);
94  const float bpos = 0.5f * ((hit0_y + hit1_y) - (hit0_x - hit1_x) * quad);
95 
96  // positive points of intersection with third layer
97  float lay2_posx = 0.0f, lay2_posy = 0.0f;
98  intersectThirdLayer(apos, bpos, hit1_x, hit1_y, lay2_posx, lay2_posy);
99 #ifdef DEBUG
100  const float lay2_posphi = getPhi(lay2_posx, lay2_posy);
101 #endif
102 
103  std::vector<int> cand_hit2_indices;
104  // MIMI lay2_hits.selectHitIndices((2.0f*hit1_z-hit0_z),(lay2_posphi+lay2_negphi)/2.0f,
105  // MIMI seed_z2cut,(lay2_posphi-lay2_negphi)/2.0f,
106  // MIMI cand_hit2_indices,true,false);
107 
108  dprint(" ihit0: " << ihit0 << " mcTrackID: " << hit0.mcTrackID(ev->simHitsInfo_) << " phi: " << hit0.phi()
109  << " z: " << hit0.z());
110  dprint(" predphi: " << (lay2_posphi + lay2_negphi) / 2.0f << "+/-" << (lay2_posphi - lay2_negphi) / 2.0f
111  << " predz: " << 2.0f * hit1_z - hit0_z << "+/-" << seed_z2cut << std::endl);
112 
113  // loop over candidate third layer hits
114  //temp_thr_seed_idcs.reserve(temp_thr_seed_idcs.size()+cand_hit2_indices.size());
115 #pragma omp simd
116  for (size_t idx = 0; idx < cand_hit2_indices.size(); ++idx) {
117  const int ihit2 = cand_hit2_indices[idx];
118  const Hit& hit2 = lay2_hits.refHit(ihit2);
119 
120  const float lay1_predz = (hit0_z + hit2.z()) / 2.0f;
121  // filter by residual of second layer hit
122  if (std::abs(lay1_predz - hit1_z) > Config::seed_z1cut)
123  continue;
124 
125  const float hit2_x = hit2.x();
126  const float hit2_y = hit2.y();
127 
128  // now fit a circle, extract pT and d0 from center and radius
129  const float mr = (hit1_y - hit0_y) / (hit1_x - hit0_x);
130  const float mt = (hit2_y - hit1_y) / (hit2_x - hit1_x);
131  const float a = (mr * mt * (hit2_y - hit0_y) + mr * (hit1_x + hit2_x) - mt * (hit0_x + hit1_x)) /
132  (2.0f * (mr - mt));
133  const float b = -1.0f * (a - (hit0_x + hit1_x) / 2.0f) / mr + (hit0_y + hit1_y) / 2.0f;
134  const float r = getHypot(hit0_x - a, hit0_y - b);
135 
136  // filter by d0 cut 5mm, pT cut 0.5 GeV (radius of 0.5 GeV track)
138  continue;
139 
140  dprint(" ihit2: " << ihit2 << " mcTrackID: " << hit2.mcTrackID(ev->simHitsInfo_)
141  << " phi: " << hit2.phi() << " z: " << hit2.z());
142 
143  temp_thr_seed_idcs.emplace_back(TripletIdx{{ihit0, ihit1, ihit2}});
144  } // end loop over third layer matches
145  } // end loop over first layer matches
146  } // end chunk of hits for parallel for
147  seed_idcs.grow_by(temp_thr_seed_idcs.begin(), temp_thr_seed_idcs.end());
148  }); // end parallel for loop over second layer hits
149  }
constexpr float dEtaSeedTrip
double z
global z - AlignmentGeometry::z0, mm
Definition: HitCollection.h:27
constexpr float lay01angdiff
void intersectThirdLayer(const float a, const float b, const float hit1_x, const float hit1_y, float &lay2_x, float &lay2_y)
Definition: seedtestMPlex.cc:9
std::array< int, 3 > TripletIdx
Definition: TrackExtra.h:41
constexpr float seed_d0cut
T sqrt(T t)
Definition: SSEVec.h:19
float getRad2(float x, float y)
Definition: Hit.h:30
Tan< T >::type tan(const T &t)
Definition: Tan.h:22
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
double f[11][100]
float getHypot(float x, float y)
Definition: Hit.h:47
#define debug
Definition: HDRShower.cc:19
constexpr float seed_z1cut
double b
Definition: hdecay.h:120
float getPhi(float x, float y)
Definition: Hit.h:34
#define dprint(x)
Definition: Debug.h:95
double a
Definition: hdecay.h:121
constexpr float seed_z0cut
constexpr float maxCurvR
std::vector< TripletIdx > TripletIdxVec
Definition: TrackExtra.h:42

◆ from_json() [1/12]

void mkfit::from_json ( const nlohmann::json nlohmann_json_j,
mkfit::LayerControl nlohmann_json_t 
)
inline

◆ from_json() [2/12]

void mkfit::from_json ( const nlohmann::ordered_json &  nlohmann_json_j,
mkfit::LayerControl nlohmann_json_t 
)
inline

Definition at line 37 of file IterationConfig.cc.

◆ from_json() [3/12]

void mkfit::from_json ( const nlohmann::json nlohmann_json_j,
mkfit::SteeringParams nlohmann_json_t 
)
inline

Definition at line 45 of file IterationConfig.cc.

◆ from_json() [4/12]

void mkfit::from_json ( const nlohmann::ordered_json &  nlohmann_json_j,
mkfit::SteeringParams nlohmann_json_t 
)
inline

Definition at line 45 of file IterationConfig.cc.

◆ from_json() [5/12]

void mkfit::from_json ( const nlohmann::json nlohmann_json_j,
mkfit::IterationLayerConfig nlohmann_json_t 
)
inline

Definition at line 54 of file IterationConfig.cc.

◆ from_json() [6/12]

void mkfit::from_json ( const nlohmann::ordered_json &  nlohmann_json_j,
mkfit::IterationLayerConfig nlohmann_json_t 
)
inline

Definition at line 54 of file IterationConfig.cc.

◆ from_json() [7/12]

void mkfit::from_json ( const nlohmann::json nlohmann_json_j,
mkfit::IterationParams nlohmann_json_t 
)
inline

Definition at line 68 of file IterationConfig.cc.

◆ from_json() [8/12]

void mkfit::from_json ( const nlohmann::ordered_json &  nlohmann_json_j,
mkfit::IterationParams nlohmann_json_t 
)
inline

Definition at line 68 of file IterationConfig.cc.

◆ from_json() [9/12]

void mkfit::from_json ( const nlohmann::json nlohmann_json_j,
mkfit::IterationConfig nlohmann_json_t 
)
inline

Definition at line 101 of file IterationConfig.cc.

◆ from_json() [10/12]

void mkfit::from_json ( const nlohmann::ordered_json &  nlohmann_json_j,
mkfit::IterationConfig nlohmann_json_t 
)
inline

Definition at line 101 of file IterationConfig.cc.

◆ from_json() [11/12]

void mkfit::from_json ( const nlohmann::json nlohmann_json_j,
mkfit::IterationsInfo nlohmann_json_t 
)
inline

Definition at line 104 of file IterationConfig.cc.

110 {

◆ from_json() [12/12]

void mkfit::from_json ( const nlohmann::ordered_json &  nlohmann_json_j,
mkfit::IterationsInfo nlohmann_json_t 
)
inline

Definition at line 104 of file IterationConfig.cc.

110 {

◆ getEta() [1/3]

float mkfit::getEta ( float  r,
float  z 
)
inline

◆ getEta() [2/3]

float mkfit::getEta ( float  theta)
inline

Definition at line 40 of file Hit.h.

References f, dqm-mbProfile::log, funct::tan(), and theta().

40 { return -1.0f * std::log(std::tan(theta / 2.0f)); }
Tan< T >::type tan(const T &t)
Definition: Tan.h:22
double f[11][100]
Geom::Theta< T > theta() const

◆ getEta() [3/3]

float mkfit::getEta ( float  x,
float  y,
float  z 
)
inline

Definition at line 42 of file Hit.h.

References f, dqm-mbProfile::log, mathSSE::sqrt(), funct::tan(), theta(), and x.

42  {
43  const float theta = std::atan2(std::sqrt(x * x + y * y), z);
44  return -1.0f * std::log(std::tan(theta / 2.0f));
45  }
float float float z
T sqrt(T t)
Definition: SSEVec.h:19
Tan< T >::type tan(const T &t)
Definition: Tan.h:22
double f[11][100]
float x
Geom::Theta< T > theta() const

◆ getEtaErr2()

float mkfit::getEtaErr2 ( float  x,
float  y,
float  z,
float  exx,
float  eyy,
float  ezz,
float  exy,
float  exz,
float  eyz 
)
inline

Definition at line 74 of file Hit.h.

References getRad2(), mathSSE::sqrt(), and x.

Referenced by mkfit::Hit::eeta(), and mkfit::TrackState::eposEta().

74  {
75  const float rad2 = getRad2(x, y);
76  const float detadx = -x / (rad2 * std::sqrt(1 + rad2 / (z * z)));
77  const float detady = -y / (rad2 * std::sqrt(1 + rad2 / (z * z)));
78  const float detadz = 1.0f / (z * std::sqrt(1 + rad2 / (z * z)));
79  return detadx * detadx * exx + detady * detady * eyy + detadz * detadz * ezz + 2.0f * detadx * detady * exy +
80  2.0f * detadx * detadz * exz + 2.0f * detady * detadz * eyz;
81  }
float float float z
T sqrt(T t)
Definition: SSEVec.h:19
float getRad2(float x, float y)
Definition: Hit.h:30
float x

◆ getHypot()

float mkfit::getHypot ( float  x,
float  y 
)
inline

Definition at line 47 of file Hit.h.

References mathSSE::sqrt(), and x.

Referenced by findSeedsByRoadSearch(), intersectThirdLayer(), mkfit::TrackState::posR(), mkfit::TrackBase::posR(), and mkfit::TrackBase::swimPhiToR().

47 { return std::sqrt(x * x + y * y); }
T sqrt(T t)
Definition: SSEVec.h:19
float x

◆ getInvRad2()

float mkfit::getInvRad2 ( float  x,
float  y 
)
inline

Definition at line 32 of file Hit.h.

References x.

32 { return 1.0f / (x * x + y * y); }
float x

◆ getInvRadErr2()

float mkfit::getInvRadErr2 ( float  x,
float  y,
float  exx,
float  eyy,
float  exy 
)
inline

Definition at line 53 of file Hit.h.

References cube(), f, getRad2(), and x.

53  {
54  return (x * x * exx + y * y * eyy + 2.0f * x * y * exy) / cube(getRad2(x, y));
55  }
float getRad2(float x, float y)
Definition: Hit.h:30
double f[11][100]
T cube(T x)
Definition: Hit.h:18
float x

◆ getMatchBin()

int mkfit::getMatchBin ( const float  pt)
inline

Definition at line 237 of file TrackExtra.cc.

References f, and DiDispStaMuonMonitor_cfi::pt.

Referenced by mkfit::TrackExtra::setCMSSWTrackIDInfoByTrkParams().

237  {
238  if (pt < 0.75f)
239  return 0;
240  else if (pt < 1.f)
241  return 1;
242  else if (pt < 2.f)
243  return 2;
244  else if (pt < 5.f)
245  return 3;
246  else if (pt < 10.f)
247  return 4;
248  else
249  return 5;
250  }
double f[11][100]

◆ getPhi()

float mkfit::getPhi ( float  x,
float  y 
)
inline

◆ getPhiErr2()

float mkfit::getPhiErr2 ( float  x,
float  y,
float  exx,
float  eyy,
float  exy 
)
inline

Definition at line 57 of file Hit.h.

References f, getRad2(), and x.

Referenced by mkfit::Hit::ephi(), and mkfit::TrackState::eposPhi().

57  {
58  const float rad2 = getRad2(x, y);
59  return (y * y * exx + x * x * eyy - 2.0f * x * y * exy) / (rad2 * rad2);
60  }
float getRad2(float x, float y)
Definition: Hit.h:30
double f[11][100]
float x

◆ getPxPxErr2()

float mkfit::getPxPxErr2 ( float  ipt,
float  phi,
float  vipt,
float  vphi 
)
inline

Definition at line 83 of file Hit.h.

References funct::cos(), and funct::sin().

Referenced by mkfit::TrackState::epxpx().

83  { // ipt = 1/pT, v = variance
84  const float iipt2 = 1.0f / (ipt * ipt); //iipt = 1/(1/pT) = pT
85  const float cosP = std::cos(phi);
86  const float sinP = std::sin(phi);
87  return iipt2 * (iipt2 * cosP * cosP * vipt + sinP * sinP * vphi);
88  }
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
Cos< T >::type cos(const T &t)
Definition: Cos.h:22

◆ getPyPyErr2()

float mkfit::getPyPyErr2 ( float  ipt,
float  phi,
float  vipt,
float  vphi 
)
inline

Definition at line 90 of file Hit.h.

References funct::cos(), and funct::sin().

Referenced by mkfit::TrackState::epypy(), and mkfit::TrackState::epzpz().

90  { // ipt = 1/pT, v = variance
91  const float iipt2 = 1.0f / (ipt * ipt); //iipt = 1/(1/pT) = pT
92  const float cosP = std::cos(phi);
93  const float sinP = std::sin(phi);
94  return iipt2 * (iipt2 * sinP * sinP * vipt + cosP * cosP * vphi);
95  }
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
Cos< T >::type cos(const T &t)
Definition: Cos.h:22

◆ getPzPzErr2()

float mkfit::getPzPzErr2 ( float  ipt,
float  theta,
float  vipt,
float  vtheta 
)
inline

Definition at line 97 of file Hit.h.

References funct::sin(), funct::tan(), and theta().

97  { // ipt = 1/pT, v = variance
98  const float iipt2 = 1.0f / (ipt * ipt); //iipt = 1/(1/pT) = pT
99  const float cotT = 1.0f / std::tan(theta);
100  const float cscT = 1.0f / std::sin(theta);
101  return iipt2 * (iipt2 * cotT * cotT * vipt + cscT * cscT * cscT * cscT * vtheta);
102  }
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
Tan< T >::type tan(const T &t)
Definition: Tan.h:22
Geom::Theta< T > theta() const

◆ getRad2()

float mkfit::getRad2 ( float  x,
float  y 
)
inline

Definition at line 30 of file Hit.h.

References x.

Referenced by conformalFitMPlex(), findSeedsByRoadSearch(), getEtaErr2(), getInvRadErr2(), getPhiErr2(), getRadErr2(), and getThetaErr2().

30 { return x * x + y * y; }
float x

◆ getRadErr2()

float mkfit::getRadErr2 ( float  x,
float  y,
float  exx,
float  eyy,
float  exy 
)
inline

Definition at line 49 of file Hit.h.

References f, getRad2(), and x.

Referenced by mkfit::TrackState::eposR().

49  {
50  return (x * x * exx + y * y * eyy + 2.0f * x * y * exy) / getRad2(x, y);
51  }
float getRad2(float x, float y)
Definition: Hit.h:30
double f[11][100]
float x

◆ getScoreCand() [1/2]

float mkfit::getScoreCand ( const track_score_func score_func,
const TrackCand cand1,
bool  penalizeTailMissHits = false,
bool  inFindCandidates = false 
)
inline

Definition at line 259 of file TrackStructures.h.

References mkfit::TrackBase::chi2(), nano_mu_local_reco_cff::chi2, mkfit::TrackCand::nFoundHits(), mkfit::TrackCand::nInsideMinusOneHits(), mkfit::TrackCand::nOverlapHits(), mkfit::TrackCand::nTailMinusOneHits(), DiDispStaMuonMonitor_cfi::pt, and mkfit::TrackBase::pT().

262  {
263  int nfoundhits = cand1.nFoundHits();
264  int noverlaphits = cand1.nOverlapHits();
265  int nmisshits = cand1.nInsideMinusOneHits();
266  int ntailmisshits = penalizeTailMissHits ? cand1.nTailMinusOneHits() : 0;
267  float pt = cand1.pT();
268  float chi2 = cand1.chi2();
269  // Do not allow for chi2<0 in score calculation
270  if (chi2 < 0)
271  chi2 = 0.f;
272  return score_func(nfoundhits, ntailmisshits, noverlaphits, nmisshits, chi2, pt, inFindCandidates);
273  }

◆ getScoreCand() [2/2]

float mkfit::getScoreCand ( const track_score_func score_func,
const Track cand1,
bool  penalizeTailMissHits = false,
bool  inFindCandidates = false 
)
inline

Definition at line 615 of file Track.h.

References mkfit::TrackBase::chi2(), nano_mu_local_reco_cff::chi2, mkfit::Track::nFoundHits(), mkfit::Track::nInsideMinusOneHits(), mkfit::Track::nOverlapHits(), mkfit::Track::nTailMinusOneHits(), DiDispStaMuonMonitor_cfi::pt, and mkfit::TrackBase::pT().

Referenced by mkfit::MkFinder::bkFitOutputTracks(), mkfit::MkBuilder::find_tracks_in_layers(), mkfit::MkFinder::findCandidates(), mkfit::CombCandidate::importSeed(), mkfit::CombCandidate::mergeCandsAndBestShortOne(), mkfit::Event::print_tracks(), mkfit::StdSeq::score_tracks(), and mkfit::TTreeValidation::setTrackScoresDumbCMSSW().

618  {
619  int nfoundhits = cand1.nFoundHits();
620  int noverlaphits = cand1.nOverlapHits();
621  int nmisshits = cand1.nInsideMinusOneHits();
622  float ntailmisshits = penalizeTailMissHits ? cand1.nTailMinusOneHits() : 0;
623  float pt = cand1.pT();
624  float chi2 = cand1.chi2();
625  // Do not allow for chi2<0 in score calculation
626  if (chi2 < 0)
627  chi2 = 0.f;
628  return score_func(nfoundhits, ntailmisshits, noverlaphits, nmisshits, chi2, pt, inFindCandidates);
629  }

◆ getScoreStruct()

float mkfit::getScoreStruct ( const track_score_func score_func,
const IdxChi2List cand1 
)
inline

Definition at line 631 of file Track.h.

References mkfit::IdxChi2List::chi2, nano_mu_local_reco_cff::chi2, mkfit::IdxChi2List::nhits, mkfit::IdxChi2List::nholes, mkfit::IdxChi2List::noverlaps, mkfit::IdxChi2List::ntailholes, DiDispStaMuonMonitor_cfi::pt, and mkfit::IdxChi2List::pt.

Referenced by mkfit::MkFinder::findCandidatesCloneEngine().

631  {
632  int nfoundhits = cand1.nhits;
633  int ntailholes = cand1.ntailholes;
634  int noverlaphits = cand1.noverlaps;
635  int nmisshits = cand1.nholes;
636  float pt = cand1.pt;
637  float chi2 = cand1.chi2;
638  // Do not allow for chi2<0 in score calculation
639  if (chi2 < 0)
640  chi2 = 0.f;
641  return score_func(nfoundhits, ntailholes, noverlaphits, nmisshits, chi2, pt, true /*inFindCandidates*/);
642  }

◆ getScoreWorstPossible()

float mkfit::getScoreWorstPossible ( )
inline

Definition at line 611 of file Track.h.

References SiStripPI::max.

Referenced by mkfit::CombCandidate::reset(), and mkfit::TrackCand::resetShortTrack().

611  {
612  return -std::numeric_limits<float>::max(); // used for handling of best short track during finding
613  }

◆ getTheta()

float mkfit::getTheta ( float  r,
float  z 
)
inline

◆ getThetaErr2()

float mkfit::getThetaErr2 ( float  x,
float  y,
float  z,
float  exx,
float  eyy,
float  ezz,
float  exy,
float  exz,
float  eyz 
)
inline

Definition at line 62 of file Hit.h.

References getRad2(), mathSSE::sqrt(), and x.

63  {
64  const float rad2 = getRad2(x, y);
65  const float rad = std::sqrt(rad2);
66  const float hypot2 = rad2 + z * z;
67  const float dthetadx = x * z / (rad * hypot2);
68  const float dthetady = y * z / (rad * hypot2);
69  const float dthetadz = -rad / hypot2;
70  return dthetadx * dthetadx * exx + dthetady * dthetady * eyy + dthetadz * dthetadz * ezz +
71  2.0f * dthetadx * dthetady * exy + 2.0f * dthetadx * dthetadz * exz + 2.0f * dthetady * dthetadz * eyz;
72  }
float float float z
T sqrt(T t)
Definition: SSEVec.h:19
float getRad2(float x, float y)
Definition: Hit.h:30
float x

◆ helixAtPlane()

void mkfit::helixAtPlane ( const MPlexLV inPar,
const MPlexQI inChg,
const MPlexHV plPnt,
const MPlexHV plNrm,
MPlexQF pathL,
MPlexLV outPar,
MPlexLL errorProp,
MPlexQI outFailFlag,
const int  N_proc,
const PropagationFlags pflags 
)

Definition at line 1066 of file PropagationMPlex.cc.

References f, and NN.

Referenced by propagateHelixToPlaneMPlex().

1075  {
1076  errorProp.setVal(0.f);
1077  outFailFlag.setVal(0.f);
1078 
1079  helixAtPlane_impl(inPar, inChg, plPnt, plNrm, pathL, outPar, errorProp, outFailFlag, 0, NN, N_proc, pflags);
1080  }
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
double f[11][100]

◆ helixAtRFromIterativeCCS()

void mkfit::helixAtRFromIterativeCCS ( const MPlexLV inPar,
const MPlexQI inChg,
const MPlexQF msRad,
MPlexLV outPar,
MPlexLL errorProp,
MPlexQI outFailFlag,
const int  N_proc,
const PropagationFlags pflags 
)

Definition at line 480 of file PropagationMPlex.cc.

References f, and NN.

Referenced by propagateHelixToRMPlex().

487  {
488  errorProp.setVal(0.f);
489  outFailFlag.setVal(0.f);
490 
491  helixAtRFromIterativeCCS_impl(inPar, inChg, msRad, outPar, errorProp, outFailFlag, 0, NN, N_proc, pflags);
492  }
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
double f[11][100]

◆ helixAtRFromIterativeCCSFullJac() [1/2]

void mkfit::helixAtRFromIterativeCCSFullJac ( const MPlexLV inPar,
const MPlexQI inChg,
const MPlexQF msRad,
MPlexLV outPar,
MPlexLL errorProp,
MPlexQI outFailFlag,
const int  N_proc 
)

◆ helixAtRFromIterativeCCSFullJac() [2/2]

void mkfit::helixAtRFromIterativeCCSFullJac ( const MPlexLV inPar,
const MPlexQI inChg,
const MPlexQF msRad,
MPlexLV outPar,
MPlexLL errorProp,
const int  N_proc 
)

Definition at line 292 of file PropagationMPlex.cc.

References funct::abs(), mkfit::Config::Bfield, CMS_UNROLL_LOOP_COUNT, funct::cos(), gather_cfg::cout, debug, dprint_np, f, g_debug, hipo(), mps_fire::i, dqmdumpme::k, dqmiodumpmetadata::n, mkfit::Config::Niter, NN, funct::sin(), sincos4(), mkfit::Const::sol, mathSSE::sqrt(), funct::tan(), theta(), and mkfit::Config::useTrigApprox.

297  {
298  errorProp.setVal(0.f);
299  MPlexLL errorPropTmp(0.f); //initialize to zero
300  MPlexLL errorPropSwap(0.f); //initialize to zero
301 
302  // loop does not vectorize with llvm16, and it issues a warning
303  // that apparently can't be suppressed with a pragma. Needs to
304  // be rechecked if future llvm versions improve vectorization.
305 #if !defined(__clang__)
306 #pragma omp simd
307 #endif
308  for (int n = 0; n < NN; ++n) {
309  //initialize erroProp to identity matrix
310  errorProp(n, 0, 0) = 1.f;
311  errorProp(n, 1, 1) = 1.f;
312  errorProp(n, 2, 2) = 1.f;
313  errorProp(n, 3, 3) = 1.f;
314  errorProp(n, 4, 4) = 1.f;
315  errorProp(n, 5, 5) = 1.f;
316 
317  const float k = inChg.constAt(n, 0, 0) * 100.f / (-Const::sol * Config::Bfield);
318  const float r = msRad.constAt(n, 0, 0);
319  float r0 = hipo(inPar.constAt(n, 0, 0), inPar.constAt(n, 1, 0));
320 
321  if (std::abs(r - r0) < 0.0001f) {
322  dprint_np(n, "distance less than 1mum, skip");
323  continue;
324  }
325 
326  const float ipt = inPar.constAt(n, 3, 0);
327  const float phiin = inPar.constAt(n, 4, 0);
328  const float theta = inPar.constAt(n, 5, 0);
329 
330  //set those that are 1. before iterations
331  errorPropTmp(n, 2, 2) = 1.f;
332  errorPropTmp(n, 3, 3) = 1.f;
333  errorPropTmp(n, 4, 4) = 1.f;
334  errorPropTmp(n, 5, 5) = 1.f;
335 
336  float cosah = 0., sinah = 0.;
337  //no trig approx here, phi and theta can be large
338  float cosP = std::cos(phiin), sinP = std::sin(phiin);
339  const float cosT = std::cos(theta), sinT = std::sin(theta);
340  float pxin = cosP / ipt;
341  float pyin = sinP / ipt;
342 
344  for (int i = 0; i < Config::Niter; ++i) {
345  dprint_np(n,
346  std::endl
347  << "attempt propagation from r=" << r0 << " to r=" << r << std::endl
348  << "x=" << outPar.At(n, 0, 0) << " y=" << outPar.At(n, 1, 0) << " z=" << outPar.At(n, 2, 0)
349  << " px=" << std::cos(phiin) / ipt << " py=" << std::sin(phiin) / ipt
350  << " pz=" << 1.f / (ipt * tan(theta)) << " q=" << inChg.constAt(n, 0, 0) << std::endl);
351 
352  r0 = hipo(outPar.constAt(n, 0, 0), outPar.constAt(n, 1, 0));
353  const float ialpha = (r - r0) * ipt / k;
354  //alpha+=ialpha;
355 
356  if constexpr (Config::useTrigApprox) {
357  sincos4(ialpha * 0.5f, sinah, cosah);
358  } else {
359  cosah = std::cos(ialpha * 0.5f);
360  sinah = std::sin(ialpha * 0.5f);
361  }
362  const float cosa = 1.f - 2.f * sinah * sinah;
363  const float sina = 2.f * sinah * cosah;
364 
365  //derivatives of alpha
366  const float dadx = -outPar.At(n, 0, 0) * ipt / (k * r0);
367  const float dady = -outPar.At(n, 1, 0) * ipt / (k * r0);
368  const float dadipt = (r - r0) / k;
369 
370  outPar.At(n, 0, 0) = outPar.constAt(n, 0, 0) + 2.f * k * sinah * (pxin * cosah - pyin * sinah);
371  outPar.At(n, 1, 0) = outPar.constAt(n, 1, 0) + 2.f * k * sinah * (pyin * cosah + pxin * sinah);
372  const float pxinold = pxin; //copy before overwriting
373  pxin = pxin * cosa - pyin * sina;
374  pyin = pyin * cosa + pxinold * sina;
375 
376  //need phi at origin, so this goes before redefining phi
377  //no trig approx here, phi can be large
378  cosP = std::cos(outPar.At(n, 4, 0));
379  sinP = std::sin(outPar.At(n, 4, 0));
380 
381  outPar.At(n, 2, 0) = outPar.constAt(n, 2, 0) + k * ialpha * cosT / (ipt * sinT);
382  outPar.At(n, 3, 0) = ipt;
383  outPar.At(n, 4, 0) = outPar.constAt(n, 4, 0) + ialpha;
384  outPar.At(n, 5, 0) = theta;
385 
386  errorPropTmp(n, 0, 0) = 1.f + k * (cosP * dadx * cosa - sinP * dadx * sina) / ipt;
387  errorPropTmp(n, 0, 1) = k * (cosP * dady * cosa - sinP * dady * sina) / ipt;
388  errorPropTmp(n, 0, 3) =
389  k * (cosP * (ipt * dadipt * cosa - sina) + sinP * ((1.f - cosa) - ipt * dadipt * sina)) / (ipt * ipt);
390  errorPropTmp(n, 0, 4) = -k * (sinP * sina + cosP * (1.f - cosa)) / ipt;
391 
392  errorPropTmp(n, 1, 0) = k * (sinP * dadx * cosa + cosP * dadx * sina) / ipt;
393  errorPropTmp(n, 1, 1) = 1.f + k * (sinP * dady * cosa + cosP * dady * sina) / ipt;
394  errorPropTmp(n, 1, 3) =
395  k * (sinP * (ipt * dadipt * cosa - sina) + cosP * (ipt * dadipt * sina - (1.f - cosa))) / (ipt * ipt);
396  errorPropTmp(n, 1, 4) = k * (cosP * sina - sinP * (1.f - cosa)) / ipt;
397 
398  errorPropTmp(n, 2, 0) = k * cosT * dadx / (ipt * sinT);
399  errorPropTmp(n, 2, 1) = k * cosT * dady / (ipt * sinT);
400  errorPropTmp(n, 2, 3) = k * cosT * (ipt * dadipt - ialpha) / (ipt * ipt * sinT);
401  errorPropTmp(n, 2, 5) = -k * ialpha / (ipt * sinT * sinT);
402 
403  errorPropTmp(n, 4, 0) = dadx;
404  errorPropTmp(n, 4, 1) = dady;
405  errorPropTmp(n, 4, 3) = dadipt;
406 
407  MultHelixPropTemp(errorProp, errorPropTmp, errorPropSwap, n);
408  errorProp = errorPropSwap;
409  }
410 
411  dprint_np(
412  n,
413  "propagation end, dump parameters"
414  << std::endl
415  << "pos = " << outPar.At(n, 0, 0) << " " << outPar.At(n, 1, 0) << " " << outPar.At(n, 2, 0) << std::endl
416  << "mom = " << std::cos(outPar.At(n, 4, 0)) / outPar.At(n, 3, 0) << " "
417  << std::sin(outPar.At(n, 4, 0)) / outPar.At(n, 3, 0) << " "
418  << 1. / (outPar.At(n, 3, 0) * tan(outPar.At(n, 5, 0)))
419  << " r=" << std::sqrt(outPar.At(n, 0, 0) * outPar.At(n, 0, 0) + outPar.At(n, 1, 0) * outPar.At(n, 1, 0))
420  << " pT=" << 1. / std::abs(outPar.At(n, 3, 0)) << std::endl);
421 
422 #ifdef DEBUG
423  if (debug && g_debug && n < N_proc) {
424  dmutex_guard;
425  std::cout << n << " jacobian" << std::endl;
426  printf("%5f %5f %5f %5f %5f %5f\n",
427  errorProp(n, 0, 0),
428  errorProp(n, 0, 1),
429  errorProp(n, 0, 2),
430  errorProp(n, 0, 3),
431  errorProp(n, 0, 4),
432  errorProp(n, 0, 5));
433  printf("%5f %5f %5f %5f %5f %5f\n",
434  errorProp(n, 1, 0),
435  errorProp(n, 1, 1),
436  errorProp(n, 1, 2),
437  errorProp(n, 1, 3),
438  errorProp(n, 1, 4),
439  errorProp(n, 1, 5));
440  printf("%5f %5f %5f %5f %5f %5f\n",
441  errorProp(n, 2, 0),
442  errorProp(n, 2, 1),
443  errorProp(n, 2, 2),
444  errorProp(n, 2, 3),
445  errorProp(n, 2, 4),
446  errorProp(n, 2, 5));
447  printf("%5f %5f %5f %5f %5f %5f\n",
448  errorProp(n, 3, 0),
449  errorProp(n, 3, 1),
450  errorProp(n, 3, 2),
451  errorProp(n, 3, 3),
452  errorProp(n, 3, 4),
453  errorProp(n, 3, 5));
454  printf("%5f %5f %5f %5f %5f %5f\n",
455  errorProp(n, 4, 0),
456  errorProp(n, 4, 1),
457  errorProp(n, 4, 2),
458  errorProp(n, 4, 3),
459  errorProp(n, 4, 4),
460  errorProp(n, 4, 5));
461  printf("%5f %5f %5f %5f %5f %5f\n",
462  errorProp(n, 5, 0),
463  errorProp(n, 5, 1),
464  errorProp(n, 5, 2),
465  errorProp(n, 5, 3),
466  errorProp(n, 5, 4),
467  errorProp(n, 5, 5));
468  }
469 #endif
470  }
471  }
Matriplex::Matriplex< float, LL, LL, NN > MPlexLL
Definition: Matrix.h:48
#define dprint_np(n, x)
Definition: Debug.h:96
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
#define CMS_UNROLL_LOOP_COUNT(N)
Definition: CMSUnrollLoop.h:48
bool g_debug
Definition: Debug.cc:2
T sqrt(T t)
Definition: SSEVec.h:19
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
constexpr float Bfield
Definition: Config.h:55
Tan< T >::type tan(const T &t)
Definition: Tan.h:22
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
double f[11][100]
constexpr bool useTrigApprox
Definition: Config.h:50
#define debug
Definition: HDRShower.cc:19
constexpr float sol
Definition: Config.h:13
float hipo(float x, float y)
Definition: Matrix.h:9
void sincos4(const float x, float &sin, float &cos)
Definition: Matrix.h:13
Geom::Theta< T > theta() const
constexpr int Niter
Definition: Config.h:49

◆ helixAtZ()

void mkfit::helixAtZ ( const MPlexLV inPar,
const MPlexQI inChg,
const MPlexQF msZ,
MPlexLV outPar,
MPlexLL errorProp,
MPlexQI outFailFlag,
const int  N_proc,
const PropagationFlags pflags 
)

Definition at line 744 of file PropagationMPlex.cc.

References funct::abs(), simBeamSpotPI::alpha, mkfit::Config::Bfield, mkfit::Config::bFieldFromZR(), funct::cos(), gather_cfg::cout, debug, l1tTrackerHTMiss_cfi::deltaZ, dprint_np, f, g_debug, hipo(), dqmdumpme::k, dqmiodumpmetadata::n, NN, DiDispStaMuonMonitor_cfi::pt, funct::sin(), sincos4(), mkfit::Const::sol, mathSSE::sqrt(), funct::tan(), theta(), mkfit::PropagationFlags::use_param_b_field, and mkfit::Config::useTrigApprox.

Referenced by propagateHelixToZMPlex().

751  {
752  errorProp.setVal(0.f);
753  outFailFlag.setVal(0.f);
754 
755  // debug = true;
756 #pragma omp simd
757  for (int n = 0; n < NN; ++n) {
758  //initialize erroProp to identity matrix, except element 2,2 which is zero
759  errorProp(n, 0, 0) = 1.f;
760  errorProp(n, 1, 1) = 1.f;
761  errorProp(n, 3, 3) = 1.f;
762  errorProp(n, 4, 4) = 1.f;
763  errorProp(n, 5, 5) = 1.f;
764  }
765  float zout[NN];
766  float zin[NN];
767  float ipt[NN];
768  float phiin[NN];
769  float theta[NN];
770 #pragma omp simd
771  for (int n = 0; n < NN; ++n) {
772  //initialize erroProp to identity matrix, except element 2,2 which is zero
773  zout[n] = msZ.constAt(n, 0, 0);
774  zin[n] = inPar.constAt(n, 2, 0);
775  ipt[n] = inPar.constAt(n, 3, 0);
776  phiin[n] = inPar.constAt(n, 4, 0);
777  theta[n] = inPar.constAt(n, 5, 0);
778  }
779 
780  float k[NN];
781  if (pflags.use_param_b_field) {
782 #pragma omp simd
783  for (int n = 0; n < NN; ++n) {
784  k[n] = inChg.constAt(n, 0, 0) * 100.f /
785  (-Const::sol * Config::bFieldFromZR(zin[n], hipo(inPar.constAt(n, 0, 0), inPar.constAt(n, 1, 0))));
786  }
787  } else {
788 #pragma omp simd
789  for (int n = 0; n < NN; ++n) {
790  k[n] = inChg.constAt(n, 0, 0) * 100.f / (-Const::sol * Config::Bfield);
791  }
792  }
793 
794  float kinv[NN];
795 #pragma omp simd
796  for (int n = 0; n < NN; ++n) {
797  kinv[n] = 1.f / k[n];
798  }
799 
800 #pragma omp simd
801  for (int n = 0; n < NN; ++n) {
802  dprint_np(n,
803  std::endl
804  << "input parameters"
805  << " inPar.constAt(n, 0, 0)=" << std::setprecision(9) << inPar.constAt(n, 0, 0)
806  << " inPar.constAt(n, 1, 0)=" << std::setprecision(9) << inPar.constAt(n, 1, 0)
807  << " inPar.constAt(n, 2, 0)=" << std::setprecision(9) << inPar.constAt(n, 2, 0)
808  << " inPar.constAt(n, 3, 0)=" << std::setprecision(9) << inPar.constAt(n, 3, 0)
809  << " inPar.constAt(n, 4, 0)=" << std::setprecision(9) << inPar.constAt(n, 4, 0)
810  << " inPar.constAt(n, 5, 0)=" << std::setprecision(9) << inPar.constAt(n, 5, 0)
811  << " inChg.constAt(n, 0, 0)=" << std::setprecision(9) << inChg.constAt(n, 0, 0));
812  }
813 #pragma omp simd
814  for (int n = 0; n < NN; ++n) {
815  dprint_np(n,
816  "propagation start, dump parameters"
817  << std::endl
818  << "pos = " << inPar.constAt(n, 0, 0) << " " << inPar.constAt(n, 1, 0) << " "
819  << inPar.constAt(n, 2, 0) << std::endl
820  << "mom (cart) = " << std::cos(inPar.constAt(n, 4, 0)) / inPar.constAt(n, 3, 0) << " "
821  << std::sin(inPar.constAt(n, 4, 0)) / inPar.constAt(n, 3, 0) << " "
822  << 1. / (inPar.constAt(n, 3, 0) * tan(inPar.constAt(n, 5, 0))) << " r="
823  << std::sqrt(inPar.constAt(n, 0, 0) * inPar.constAt(n, 0, 0) +
824  inPar.constAt(n, 1, 0) * inPar.constAt(n, 1, 0))
825  << " pT=" << 1. / std::abs(inPar.constAt(n, 3, 0)) << " q=" << inChg.constAt(n, 0, 0)
826  << " targetZ=" << msZ.constAt(n, 0, 0) << std::endl);
827  }
828 
829  float pt[NN];
830 #pragma omp simd
831  for (int n = 0; n < NN; ++n) {
832  pt[n] = 1.f / ipt[n];
833  }
834 
835  //no trig approx here, phi can be large
836  float cosP[NN];
837  float sinP[NN];
838 #pragma omp simd
839  for (int n = 0; n < NN; ++n) {
840  cosP[n] = std::cos(phiin[n]);
841  }
842 
843 #pragma omp simd
844  for (int n = 0; n < NN; ++n) {
845  sinP[n] = std::sin(phiin[n]);
846  }
847 
848  float cosT[NN];
849  float sinT[NN];
850 #pragma omp simd
851  for (int n = 0; n < NN; ++n) {
852  cosT[n] = std::cos(theta[n]);
853  }
854 
855 #pragma omp simd
856  for (int n = 0; n < NN; ++n) {
857  sinT[n] = std::sin(theta[n]);
858  }
859 
860  float tanT[NN];
861  float icos2T[NN];
862  float pxin[NN];
863  float pyin[NN];
864 #pragma omp simd
865  for (int n = 0; n < NN; ++n) {
866  tanT[n] = sinT[n] / cosT[n];
867  icos2T[n] = 1.f / (cosT[n] * cosT[n]);
868  pxin[n] = cosP[n] * pt[n];
869  pyin[n] = sinP[n] * pt[n];
870  }
871 
872  float deltaZ[NN];
873  float alpha[NN];
874 #pragma omp simd
875  for (int n = 0; n < NN; ++n) {
876  deltaZ[n] = zout[n] - zin[n];
877  alpha[n] = deltaZ[n] * tanT[n] * ipt[n] * kinv[n];
878  }
879 
880  float cosahTmp[NN];
881  float sinahTmp[NN];
882  if constexpr (Config::useTrigApprox) {
883 #if !defined(__INTEL_COMPILER)
884 #pragma omp simd
885 #endif
886  for (int n = 0; n < NN; ++n) {
887  sincos4(alpha[n] * 0.5f, sinahTmp[n], cosahTmp[n]);
888  }
889  } else {
890 #if !defined(__INTEL_COMPILER)
891 #pragma omp simd
892 #endif
893  for (int n = 0; n < NN; ++n) {
894  cosahTmp[n] = std::cos(alpha[n] * 0.5f);
895  }
896 #if !defined(__INTEL_COMPILER)
897 #pragma omp simd
898 #endif
899  for (int n = 0; n < NN; ++n) {
900  sinahTmp[n] = std::sin(alpha[n] * 0.5f);
901  }
902  }
903 
904  float cosah[NN];
905  float sinah[NN];
906  float cosa[NN];
907  float sina[NN];
908 #pragma omp simd
909  for (int n = 0; n < NN; ++n) {
910  cosah[n] = cosahTmp[n];
911  sinah[n] = sinahTmp[n];
912  cosa[n] = 1.f - 2.f * sinah[n] * sinah[n];
913  sina[n] = 2.f * sinah[n] * cosah[n];
914  }
915 
916 //update parameters
917 #pragma omp simd
918  for (int n = 0; n < NN; ++n) {
919  outPar.At(n, 0, 0) = outPar.At(n, 0, 0) + 2.f * k[n] * sinah[n] * (pxin[n] * cosah[n] - pyin[n] * sinah[n]);
920  outPar.At(n, 1, 0) = outPar.At(n, 1, 0) + 2.f * k[n] * sinah[n] * (pyin[n] * cosah[n] + pxin[n] * sinah[n]);
921  outPar.At(n, 2, 0) = zout[n];
922  outPar.At(n, 4, 0) = phiin[n] + alpha[n];
923  }
924 
925 #pragma omp simd
926  for (int n = 0; n < NN; ++n) {
927  dprint_np(n,
928  "propagation to Z end (OLD), dump parameters\n"
929  << " pos = " << outPar(n, 0, 0) << " " << outPar(n, 1, 0) << " " << outPar(n, 2, 0) << "\t\t r="
930  << std::sqrt(outPar(n, 0, 0) * outPar(n, 0, 0) + outPar(n, 1, 0) * outPar(n, 1, 0)) << std::endl
931  << " mom = " << outPar(n, 3, 0) << " " << outPar(n, 4, 0) << " " << outPar(n, 5, 0) << std::endl
932  << " cart= " << std::cos(outPar(n, 4, 0)) / outPar(n, 3, 0) << " "
933  << std::sin(outPar(n, 4, 0)) / outPar(n, 3, 0) << " "
934  << 1. / (outPar(n, 3, 0) * tan(outPar(n, 5, 0))) << "\t\tpT=" << 1. / std::abs(outPar(n, 3, 0))
935  << std::endl);
936  }
937 
938  float pxcaMpysa[NN];
939 #pragma omp simd
940  for (int n = 0; n < NN; ++n) {
941  pxcaMpysa[n] = pxin[n] * cosa[n] - pyin[n] * sina[n];
942  }
943 
944 #pragma omp simd
945  for (int n = 0; n < NN; ++n) {
946  errorProp(n, 0, 2) = -tanT[n] * ipt[n] * pxcaMpysa[n];
947  errorProp(n, 0, 3) =
948  k[n] * pt[n] * pt[n] *
949  (cosP[n] * (alpha[n] * cosa[n] - sina[n]) + sinP[n] * 2.f * sinah[n] * (sinah[n] - alpha[n] * cosah[n]));
950  errorProp(n, 0, 4) = -2.f * k[n] * pt[n] * sinah[n] * (sinP[n] * cosah[n] + cosP[n] * sinah[n]);
951  errorProp(n, 0, 5) = deltaZ[n] * ipt[n] * pxcaMpysa[n] * icos2T[n];
952  }
953 
954  float pycaPpxsa[NN];
955 #pragma omp simd
956  for (int n = 0; n < NN; ++n) {
957  pycaPpxsa[n] = pyin[n] * cosa[n] + pxin[n] * sina[n];
958  }
959 
960 #pragma omp simd
961  for (int n = 0; n < NN; ++n) {
962  errorProp(n, 1, 2) = -tanT[n] * ipt[n] * pycaPpxsa[n];
963  errorProp(n, 1, 3) =
964  k[n] * pt[n] * pt[n] *
965  (sinP[n] * (alpha[n] * cosa[n] - sina[n]) - cosP[n] * 2.f * sinah[n] * (sinah[n] - alpha[n] * cosah[n]));
966  errorProp(n, 1, 4) = 2.f * k[n] * pt[n] * sinah[n] * (cosP[n] * cosah[n] - sinP[n] * sinah[n]);
967  errorProp(n, 1, 5) = deltaZ[n] * ipt[n] * pycaPpxsa[n] * icos2T[n];
968  }
969 
970 #pragma omp simd
971  for (int n = 0; n < NN; ++n) {
972  errorProp(n, 4, 2) = -ipt[n] * tanT[n] * kinv[n];
973  errorProp(n, 4, 3) = tanT[n] * deltaZ[n] * kinv[n];
974  errorProp(n, 4, 5) = ipt[n] * deltaZ[n] * kinv[n] * icos2T[n];
975  }
976 
977 #pragma omp simd
978  for (int n = 0; n < NN; ++n) {
979  dprint_np(
980  n,
981  "propagation end, dump parameters"
982  << std::endl
983  << "pos = " << outPar.At(n, 0, 0) << " " << outPar.At(n, 1, 0) << " " << outPar.At(n, 2, 0) << std::endl
984  << "mom (cart) = " << std::cos(outPar.At(n, 4, 0)) / outPar.At(n, 3, 0) << " "
985  << std::sin(outPar.At(n, 4, 0)) / outPar.At(n, 3, 0) << " "
986  << 1. / (outPar.At(n, 3, 0) * tan(outPar.At(n, 5, 0)))
987  << " r=" << std::sqrt(outPar.At(n, 0, 0) * outPar.At(n, 0, 0) + outPar.At(n, 1, 0) * outPar.At(n, 1, 0))
988  << " pT=" << 1. / std::abs(outPar.At(n, 3, 0)) << std::endl);
989  }
990 
991  // PROP-FAIL-ENABLE Disabled to keep physics changes minimal.
992  // To be reviewed, enabled and processed accordingly elsewhere.
993  /*
994  // Check for errors, set fail-flag.
995  for (int n = 0; n < NN; ++n) {
996  // We propagate for alpha: mark fail when prop angle more than pi/2
997  if (std::abs(alpha[n]) > 1.57) {
998  dprintf("helixAtZ: more than quarter turn, alpha = %f\n", alpha[n]);
999  outFailFlag[n] = 1;
1000  } else {
1001  // Have we reached desired z? We can't know, we copy desired z to actual z.
1002  // Are we close to apex? Same condition as in propToR, 12.5 deg, cos(78.5deg) = 0.2
1003  float dotp = (outPar.At(n, 0, 0) * std::cos(outPar.At(n, 4, 0)) +
1004  outPar.At(n, 1, 0) * std::sin(outPar.At(n, 4, 0))) /
1005  std::hypot(outPar.At(n, 0, 0), outPar.At(n, 1, 0));
1006  if (dotp < 0.2 || dotp < 0) {
1007  dprintf("helixAtZ: dot product bad, dotp = %f\n", dotp);
1008  outFailFlag[n] = 1;
1009  }
1010  }
1011  }
1012  */
1013 
1014 #ifdef DEBUG
1015  if (debug && g_debug) {
1016  for (int n = 0; n < N_proc; ++n) {
1017  dmutex_guard;
1018  std::cout << n << ": jacobian" << std::endl;
1019  printf("%5f %5f %5f %5f %5f %5f\n",
1020  errorProp(n, 0, 0),
1021  errorProp(n, 0, 1),
1022  errorProp(n, 0, 2),
1023  errorProp(n, 0, 3),
1024  errorProp(n, 0, 4),
1025  errorProp(n, 0, 5));
1026  printf("%5f %5f %5f %5f %5f %5f\n",
1027  errorProp(n, 1, 0),
1028  errorProp(n, 1, 1),
1029  errorProp(n, 1, 2),
1030  errorProp(n, 1, 3),
1031  errorProp(n, 1, 4),
1032  errorProp(n, 1, 5));
1033  printf("%5f %5f %5f %5f %5f %5f\n",
1034  errorProp(n, 2, 0),
1035  errorProp(n, 2, 1),
1036  errorProp(n, 2, 2),
1037  errorProp(n, 2, 3),
1038  errorProp(n, 2, 4),
1039  errorProp(n, 2, 5));
1040  printf("%5f %5f %5f %5f %5f %5f\n",
1041  errorProp(n, 3, 0),
1042  errorProp(n, 3, 1),
1043  errorProp(n, 3, 2),
1044  errorProp(n, 3, 3),
1045  errorProp(n, 3, 4),
1046  errorProp(n, 3, 5));
1047  printf("%5f %5f %5f %5f %5f %5f\n",
1048  errorProp(n, 4, 0),
1049  errorProp(n, 4, 1),
1050  errorProp(n, 4, 2),
1051  errorProp(n, 4, 3),
1052  errorProp(n, 4, 4),
1053  errorProp(n, 4, 5));
1054  printf("%5f %5f %5f %5f %5f %5f\n",
1055  errorProp(n, 5, 0),
1056  errorProp(n, 5, 1),
1057  errorProp(n, 5, 2),
1058  errorProp(n, 5, 3),
1059  errorProp(n, 5, 4),
1060  errorProp(n, 5, 5));
1061  }
1062  }
1063 #endif
1064  }
#define dprint_np(n, x)
Definition: Debug.h:96
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
bool g_debug
Definition: Debug.cc:2
T sqrt(T t)
Definition: SSEVec.h:19
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
constexpr float Bfield
Definition: Config.h:55
Tan< T >::type tan(const T &t)
Definition: Tan.h:22
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
double f[11][100]
constexpr bool useTrigApprox
Definition: Config.h:50
#define debug
Definition: HDRShower.cc:19
constexpr float sol
Definition: Config.h:13
float hipo(float x, float y)
Definition: Matrix.h:9
float bFieldFromZR(const float z, const float r)
Definition: Config.h:126
void sincos4(const float x, float &sin, float &cos)
Definition: Matrix.h:13
Geom::Theta< T > theta() const

◆ hipo()

float mkfit::hipo ( float  x,
float  y 
)
inline

◆ hipo_sqr()

float mkfit::hipo_sqr ( float  x,
float  y 
)
inline

Definition at line 11 of file Matrix.h.

References x.

Referenced by mkfit::MkBase::radiusSqr().

11 { return x * x + y * y; }
float x

◆ intersectThirdLayer()

void mkfit::intersectThirdLayer ( const float  a,
const float  b,
const float  hit1_x,
const float  hit1_y,
float &  lay2_x,
float &  lay2_y 
)
inline

Definition at line 9 of file seedtestMPlex.cc.

References a, testProducerWithPsetDescEmpty_cfi::a2, b, b2, f, getHypot(), mkfit::Config::maxCurvR, and mathSSE::sqrt().

Referenced by findSeedsByRoadSearch().

10  {
11  const float a2 = a * a;
12  const float b2 = b * b;
13  const float a2b2 = a2 + b2;
14  const float lay2rad2 = (Config::fRadialSpacing * Config::fRadialSpacing) * 9.0f; // average third radius squared
15  const float maxCurvR2 = Config::maxCurvR * Config::maxCurvR;
16 
17  const float quad =
18  std::sqrt(2.0f * maxCurvR2 * (a2b2 + lay2rad2) - (a2b2 - lay2rad2) * (a2b2 - lay2rad2) - maxCurvR2 * maxCurvR2);
19  const float pos[2] = {(a2 * a + a * (b2 + lay2rad2 - maxCurvR2) - b * quad) / a2b2,
20  (b2 * b + b * (a2 + lay2rad2 - maxCurvR2) + a * quad) / a2b2};
21  const float neg[2] = {(a2 * a + a * (b2 + lay2rad2 - maxCurvR2) + b * quad) / a2b2,
22  (b2 * b + b * (a2 + lay2rad2 - maxCurvR2) - a * quad) / a2b2};
23 
24  // since we have two intersection points, arbitrate which one is closer to layer2 hit
25  if (getHypot(pos[0] - hit1_x, pos[1] - hit1_y) < getHypot(neg[0] - hit1_x, neg[1] - hit1_y)) {
26  lay2_x = pos[0];
27  lay2_y = pos[1];
28  } else {
29  lay2_x = neg[0];
30  lay2_y = neg[1];
31  }
32  }
T sqrt(T t)
Definition: SSEVec.h:19
double f[11][100]
float getHypot(float x, float y)
Definition: Hit.h:47
weight_default_t b2[10]
Definition: b2.h:9
double b
Definition: hdecay.h:120
double a
Definition: hdecay.h:121
constexpr float maxCurvR

◆ isFinite()

constexpr bool mkfit::isFinite ( float  x)

Definition at line 13 of file cms_common_macros.h.

References ztail::d, edm::isFinite(), MainPageGenerator::l, gpuClustering::pixelStatus::mask, findQualityFiles::v, and x.

Referenced by mkfit::MkFinder::bkFitOutputTracks(), mkfit::TrackBase::hasNanNSillyValues(), mkfit::TrackBase::hasSillyValues(), and mkfit::MkFinder::selectHitIndices().

13  {
14 #ifdef MKFIT_STANDALONE
15  const unsigned int mask = 0x7f800000;
16  union {
17  unsigned int l;
18  float d;
19  } v = {.d = x};
20  return (v.l & mask) != mask;
21 #else
22  return edm::isFinite(x);
23 #endif
24  }
constexpr bool isFinite(T x)
constexpr uint32_t mask
Definition: gpuClustering.h:26
d
Definition: ztail.py:151
float x

◆ isStripQCompatible()

bool mkfit::isStripQCompatible ( int  itrack,
bool  isBarrel,
const MPlexLS pErr,
const MPlexLV pPar,
const MPlexHS msErr,
const MPlexHV msPar 
)

Definition at line 1185 of file MkFinder.cc.

References funct::abs(), dprint, f, PixelPluginsPhase0_cfi::isBarrel, SiStripPI::max, amptDefault_cfi::proj, and mathSSE::sqrt().

Referenced by mkfit::MkFinder::findCandidates(), and mkfit::MkFinder::findCandidatesCloneEngine().

1186  {
1187  //check module compatibility via long strip side = L/sqrt(12)
1188  if (isBarrel) { //check z direction only
1189  const float res = std::abs(msPar.constAt(itrack, 2, 0) - pPar.constAt(itrack, 2, 0));
1190  const float hitHL = sqrt(msErr.constAt(itrack, 2, 2) * 3.f); //half-length
1191  const float qErr = sqrt(pErr.constAt(itrack, 2, 2));
1192  dprint("qCompat " << hitHL << " + " << 3.f * qErr << " vs " << res);
1193  return hitHL + std::max(3.f * qErr, 0.5f) > res;
1194  } else { //project on xy, assuming the strip Length >> Width
1195  const float res[2]{msPar.constAt(itrack, 0, 0) - pPar.constAt(itrack, 0, 0),
1196  msPar.constAt(itrack, 1, 0) - pPar.constAt(itrack, 1, 0)};
1197  const float hitT2 = msErr.constAt(itrack, 0, 0) + msErr.constAt(itrack, 1, 1);
1198  const float hitT2inv = 1.f / hitT2;
1199  const float proj[3] = {msErr.constAt(itrack, 0, 0) * hitT2inv,
1200  msErr.constAt(itrack, 0, 1) * hitT2inv,
1201  msErr.constAt(itrack, 1, 1) * hitT2inv};
1202  const float qErr =
1203  sqrt(std::abs(pErr.constAt(itrack, 0, 0) * proj[0] + 2.f * pErr.constAt(itrack, 0, 1) * proj[1] +
1204  pErr.constAt(itrack, 1, 1) * proj[2])); //take abs to avoid non-pos-def cases
1205  const float resProj =
1206  sqrt(res[0] * proj[0] * res[0] + 2.f * res[1] * proj[1] * res[0] + res[1] * proj[2] * res[1]);
1207  dprint("qCompat " << sqrt(hitT2 * 3.f) << " + " << 3.f * qErr << " vs " << resProj);
1208  return sqrt(hitT2 * 3.f) + std::max(3.f * qErr, 0.5f) > resProj;
1209  }
1210  }
Definition: Electron.h:6
T sqrt(T t)
Definition: SSEVec.h:19
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
double f[11][100]
#define dprint(x)
Definition: Debug.h:95

◆ kalmanComputeChi2()

void mkfit::kalmanComputeChi2 ( const MPlexLS psErr,
const MPlexLV psPar,
const MPlexQI inChg,
const MPlexHS msErr,
const MPlexHV msPar,
MPlexQF outChi2,
const int  N_proc 
)

Definition at line 753 of file KalmanUtilsMPlex.cc.

References kalmanOperation(), and KFO_Calculate_Chi2.

759  {
760  kalmanOperation(KFO_Calculate_Chi2, psErr, psPar, msErr, msPar, dummy_err, dummy_par, outChi2, N_proc);
761  }
void kalmanOperation(const int kfOp, const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, MPlexLS &outErr, MPlexLV &outPar, MPlexQF &outChi2, const int N_proc)

◆ kalmanComputeChi2Endcap()

void mkfit::kalmanComputeChi2Endcap ( const MPlexLS psErr,
const MPlexLV psPar,
const MPlexQI inChg,
const MPlexHS msErr,
const MPlexHV msPar,
MPlexQF outChi2,
const int  N_proc 
)

Definition at line 1371 of file KalmanUtilsMPlex.cc.

References kalmanOperationEndcap(), and KFO_Calculate_Chi2.

1377  {
1378  kalmanOperationEndcap(KFO_Calculate_Chi2, psErr, psPar, msErr, msPar, dummy_err, dummy_par, outChi2, N_proc);
1379  }
void kalmanOperationEndcap(const int kfOp, const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, MPlexLS &outErr, MPlexLV &outPar, MPlexQF &outChi2, const int N_proc)

◆ kalmanComputeChi2Plane()

void mkfit::kalmanComputeChi2Plane ( const MPlexLS psErr,
const MPlexLV psPar,
const MPlexQI inChg,
const MPlexHS msErr,
const MPlexHV msPar,
const MPlexHV plNrm,
const MPlexHV plDir,
MPlexQF outChi2,
const int  N_proc 
)

Definition at line 1066 of file KalmanUtilsMPlex.cc.

References kalmanOperationPlane(), and KFO_Calculate_Chi2.

1074  {
1076  KFO_Calculate_Chi2, psErr, psPar, msErr, msPar, plNrm, plDir, dummy_err, dummy_par, outChi2, N_proc);
1077  }
void kalmanOperationPlane(const int kfOp, const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, const MPlexHV &plNrm, const MPlexHV &plDir, MPlexLS &outErr, MPlexLV &outPar, MPlexQF &outChi2, const int N_proc)

◆ kalmanOperation()

void mkfit::kalmanOperation ( const int  kfOp,
const MPlexLS psErr,
const MPlexLV psPar,
const MPlexHS msErr,
const MPlexHV msPar,
MPlexLS outErr,
MPlexLV outPar,
MPlexQF outChi2,
const int  N_proc 
)

Definition at line 793 of file KalmanUtilsMPlex.cc.

References Matriplex::hypot(), mps_fire::i, Matriplex::invertCramerSym(), dqmiolumiharvest::j, KFO_Calculate_Chi2, KFO_Local_Cov, KFO_Update_Params, dqmiodumpmetadata::n, NN, and squashPhiMPlex().

Referenced by mkfit::MkFinder::bkFitFitTracks(), mkfit::MkFinder::bkFitFitTracksBH(), kalmanComputeChi2(), kalmanPropagateAndComputeChi2(), kalmanPropagateAndUpdate(), and kalmanUpdate().

801  {
802 #ifdef DEBUG
803  {
804  dmutex_guard;
805  printf("psPar:\n");
806  for (int i = 0; i < 6; ++i) {
807  printf("%8f ", psPar.constAt(0, 0, i));
808  printf("\n");
809  }
810  printf("\n");
811  printf("psErr:\n");
812  for (int i = 0; i < 6; ++i) {
813  for (int j = 0; j < 6; ++j)
814  printf("%8f ", psErr.constAt(0, i, j));
815  printf("\n");
816  }
817  printf("\n");
818  printf("msPar:\n");
819  for (int i = 0; i < 3; ++i) {
820  printf("%8f ", msPar.constAt(0, 0, i));
821  printf("\n");
822  }
823  printf("\n");
824  printf("msErr:\n");
825  for (int i = 0; i < 3; ++i) {
826  for (int j = 0; j < 3; ++j)
827  printf("%8f ", msErr.constAt(0, i, j));
828  printf("\n");
829  }
830  printf("\n");
831  }
832 #endif
833 
834  // Rotate global point on tangent plane to cylinder
835  // Tangent point is half way between hit and propagate position
836 
837  // Rotation matrix
838  // rotT00 0 rotT01
839  // rotT01 0 -rotT00
840  // 0 1 0
841  // Minimize temporaries: only two float are needed!
842 
843  MPlexQF rotT00;
844  MPlexQF rotT01;
845  for (int n = 0; n < NN; ++n) {
846  const float r = std::hypot(msPar.constAt(n, 0, 0), msPar.constAt(n, 1, 0));
847  rotT00.At(n, 0, 0) = -(msPar.constAt(n, 1, 0) + psPar.constAt(n, 1, 0)) / (2 * r);
848  rotT01.At(n, 0, 0) = (msPar.constAt(n, 0, 0) + psPar.constAt(n, 0, 0)) / (2 * r);
849  }
850 
851  MPlexHV res_glo; //position residual in global coordinates
852  SubtractFirst3(msPar, psPar, res_glo);
853 
854  MPlexHS resErr_glo; //covariance sum in global position coordinates
855  AddIntoUpperLeft3x3(psErr, msErr, resErr_glo);
856 
857  MPlex2V res_loc; //position residual in local coordinates
858  RotateResidualsOnTangentPlane(rotT00, rotT01, res_glo, res_loc);
859  MPlex2S resErr_loc; //covariance sum in local position coordinates
860  MPlexHH tempHH;
861  ProjectResErr(rotT00, rotT01, resErr_glo, tempHH);
862  ProjectResErrTransp(rotT00, rotT01, tempHH, resErr_loc);
863 
864 #ifdef DEBUG
865  {
866  dmutex_guard;
867  printf("res_glo:\n");
868  for (int i = 0; i < 3; ++i) {
869  printf("%8f ", res_glo.At(0, i, 0));
870  }
871  printf("\n");
872  printf("resErr_glo:\n");
873  for (int i = 0; i < 3; ++i) {
874  for (int j = 0; j < 3; ++j)
875  printf("%8f ", resErr_glo.At(0, i, j));
876  printf("\n");
877  }
878  printf("\n");
879  printf("res_loc:\n");
880  for (int i = 0; i < 2; ++i) {
881  printf("%8f ", res_loc.At(0, i, 0));
882  }
883  printf("\n");
884  printf("tempHH:\n");
885  for (int i = 0; i < 3; ++i) {
886  for (int j = 0; j < 3; ++j)
887  printf("%8f ", tempHH.At(0, i, j));
888  printf("\n");
889  }
890  printf("\n");
891  printf("resErr_loc:\n");
892  for (int i = 0; i < 2; ++i) {
893  for (int j = 0; j < 2; ++j)
894  printf("%8f ", resErr_loc.At(0, i, j));
895  printf("\n");
896  }
897  printf("\n");
898  }
899 #endif
900 
901  //invert the 2x2 matrix
902  Matriplex::invertCramerSym(resErr_loc);
903 
904  if (kfOp & KFO_Calculate_Chi2) {
905  Chi2Similarity(res_loc, resErr_loc, outChi2);
906 
907 #ifdef DEBUG
908  {
909  dmutex_guard;
910  printf("resErr_loc (Inv):\n");
911  for (int i = 0; i < 2; ++i) {
912  for (int j = 0; j < 2; ++j)
913  printf("%8f ", resErr_loc.At(0, i, j));
914  printf("\n");
915  }
916  printf("\n");
917  printf("chi2: %8f\n", outChi2.At(0, 0, 0));
918  }
919 #endif
920  }
921 
922  if (kfOp & KFO_Update_Params) {
923  MPlexLS psErrLoc = psErr;
924  if (kfOp & KFO_Local_Cov)
925  CovXYconstrain(rotT00, rotT01, psErr, psErrLoc);
926 
927  MPlexLH K; // kalman gain, fixme should be L2
928  KalmanHTG(rotT00, rotT01, resErr_loc, tempHH); // intermediate term to get kalman gain (H^T*G)
929  KalmanGain(psErrLoc, tempHH, K);
930 
931  MultResidualsAdd(K, psPar, res_loc, outPar);
932 
933  squashPhiMPlex(outPar, N_proc); // ensure phi is between |pi|
934 
935  MPlexLL tempLL;
936  KHMult(K, rotT00, rotT01, tempLL);
937  KHC(tempLL, psErrLoc, outErr);
938  outErr.subtract(psErrLoc, outErr);
939 
940 #ifdef DEBUG
941  {
942  dmutex_guard;
943  if (kfOp & KFO_Local_Cov) {
944  printf("psErrLoc:\n");
945  for (int i = 0; i < 6; ++i) {
946  for (int j = 0; j < 6; ++j)
947  printf("% 8e ", psErrLoc.At(0, i, j));
948  printf("\n");
949  }
950  printf("\n");
951  }
952  printf("resErr_loc (Inv):\n");
953  for (int i = 0; i < 2; ++i) {
954  for (int j = 0; j < 2; ++j)
955  printf("%8f ", resErr_loc.At(0, i, j));
956  printf("\n");
957  }
958  printf("\n");
959  printf("tempHH:\n");
960  for (int i = 0; i < 3; ++i) {
961  for (int j = 0; j < 3; ++j)
962  printf("%8f ", tempHH.At(0, i, j));
963  printf("\n");
964  }
965  printf("\n");
966  printf("K:\n");
967  for (int i = 0; i < 6; ++i) {
968  for (int j = 0; j < 3; ++j)
969  printf("%8f ", K.At(0, i, j));
970  printf("\n");
971  }
972  printf("\n");
973  printf("tempLL:\n");
974  for (int i = 0; i < 6; ++i) {
975  for (int j = 0; j < 6; ++j)
976  printf("%8f ", tempLL.At(0, i, j));
977  printf("\n");
978  }
979  printf("\n");
980  printf("outPar:\n");
981  for (int i = 0; i < 6; ++i) {
982  printf("%8f ", outPar.At(0, i, 0));
983  }
984  printf("\n");
985  printf("outErr:\n");
986  for (int i = 0; i < 6; ++i) {
987  for (int j = 0; j < 6; ++j)
988  printf("%8f ", outErr.At(0, i, j));
989  printf("\n");
990  }
991  printf("\n");
992  }
993 #endif
994  }
995  }
Matriplex::Matriplex< float, HH, HH, NN > MPlexHH
Definition: Matrix.h:52
Matriplex::Matriplex< float, LL, LL, NN > MPlexLL
Definition: Matrix.h:48
MPlex< T, D1, D2, N > hypot(const MPlex< T, D1, D2, N > &a, const MPlex< T, D1, D2, N > &b)
Definition: Matriplex.h:436
Matriplex::Matriplex< float, HH, 1, NN > MPlexHV
Definition: Matrix.h:53
void squashPhiMPlex(MPlexLV &par, const int N_proc)
void invertCramerSym(MPlexSym< T, D, N > &A, double *determ=nullptr)
Definition: MatriplexSym.h:418
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
Matriplex::Matriplex< float, 1, 1, NN > MPlexQF
Definition: Matrix.h:71
Matriplex::MatriplexSym< float, 2, NN > MPlex2S
Definition: Matrix.h:62
Matriplex::MatriplexSym< float, HH, NN > MPlexHS
Definition: Matrix.h:54
Matriplex::Matriplex< float, LL, HH, NN > MPlexLH
Definition: Matrix.h:64
Matriplex::Matriplex< float, 2, 1, NN > MPlex2V
Definition: Matrix.h:61
Matriplex::MatriplexSym< float, LL, NN > MPlexLS
Definition: Matrix.h:50

◆ kalmanOperationEndcap()

void mkfit::kalmanOperationEndcap ( const int  kfOp,
const MPlexLS psErr,
const MPlexLV psPar,
const MPlexHS msErr,
const MPlexHV msPar,
MPlexLS outErr,
MPlexLV outPar,
MPlexQF outChi2,
const int  N_proc 
)

Definition at line 1411 of file KalmanUtilsMPlex.cc.

References mps_fire::i, Matriplex::invertCramerSym(), dqmiolumiharvest::j, KFO_Calculate_Chi2, KFO_Update_Params, and squashPhiMPlex().

Referenced by mkfit::MkFinder::bkFitFitTracks(), mkfit::MkFinder::bkFitFitTracksBH(), kalmanComputeChi2Endcap(), kalmanPropagateAndComputeChi2Endcap(), kalmanPropagateAndUpdateEndcap(), and kalmanUpdateEndcap().

1419  {
1420 #ifdef DEBUG
1421  {
1422  dmutex_guard;
1423  printf("updateParametersEndcapMPlex\n");
1424  printf("psPar:\n");
1425  for (int i = 0; i < 6; ++i) {
1426  printf("%8f ", psPar.constAt(0, 0, i));
1427  printf("\n");
1428  }
1429  printf("\n");
1430  printf("msPar:\n");
1431  for (int i = 0; i < 3; ++i) {
1432  printf("%8f ", msPar.constAt(0, 0, i));
1433  printf("\n");
1434  }
1435  printf("\n");
1436  printf("psErr:\n");
1437  for (int i = 0; i < 6; ++i) {
1438  for (int j = 0; j < 6; ++j)
1439  printf("%8f ", psErr.constAt(0, i, j));
1440  printf("\n");
1441  }
1442  printf("\n");
1443  printf("msErr:\n");
1444  for (int i = 0; i < 3; ++i) {
1445  for (int j = 0; j < 3; ++j)
1446  printf("%8f ", msErr.constAt(0, i, j));
1447  printf("\n");
1448  }
1449  printf("\n");
1450  }
1451 #endif
1452 
1453  MPlex2V res;
1454  SubtractFirst2(msPar, psPar, res);
1455 
1456  MPlex2S resErr;
1457  AddIntoUpperLeft2x2(psErr, msErr, resErr);
1458 
1459 #ifdef DEBUG
1460  {
1461  dmutex_guard;
1462  printf("resErr:\n");
1463  for (int i = 0; i < 2; ++i) {
1464  for (int j = 0; j < 2; ++j)
1465  printf("%8f ", resErr.At(0, i, j));
1466  printf("\n");
1467  }
1468  printf("\n");
1469  }
1470 #endif
1471 
1472  //invert the 2x2 matrix
1474 
1475  if (kfOp & KFO_Calculate_Chi2) {
1476  Chi2Similarity(res, resErr, outChi2);
1477 
1478 #ifdef DEBUG
1479  {
1480  dmutex_guard;
1481  printf("resErr_loc (Inv):\n");
1482  for (int i = 0; i < 2; ++i) {
1483  for (int j = 0; j < 2; ++j)
1484  printf("%8f ", resErr.At(0, i, j));
1485  printf("\n");
1486  }
1487  printf("\n");
1488  printf("chi2: %8f\n", outChi2.At(0, 0, 0));
1489  }
1490 #endif
1491  }
1492 
1493  if (kfOp & KFO_Update_Params) {
1494  MPlexL2 K;
1495  KalmanGain(psErr, resErr, K);
1496 
1497  MultResidualsAdd(K, psPar, res, outPar);
1498 
1499  squashPhiMPlex(outPar, N_proc); // ensure phi is between |pi|
1500 
1501  KHC(K, psErr, outErr);
1502 
1503 #ifdef DEBUG
1504  {
1505  dmutex_guard;
1506  printf("outErr before subtract:\n");
1507  for (int i = 0; i < 6; ++i) {
1508  for (int j = 0; j < 6; ++j)
1509  printf("%8f ", outErr.At(0, i, j));
1510  printf("\n");
1511  }
1512  printf("\n");
1513  }
1514 #endif
1515 
1516  outErr.subtract(psErr, outErr);
1517 
1518 #ifdef DEBUG
1519  {
1520  dmutex_guard;
1521  printf("res:\n");
1522  for (int i = 0; i < 2; ++i) {
1523  printf("%8f ", res.At(0, i, 0));
1524  }
1525  printf("\n");
1526  printf("resErr (Inv):\n");
1527  for (int i = 0; i < 2; ++i) {
1528  for (int j = 0; j < 2; ++j)
1529  printf("%8f ", resErr.At(0, i, j));
1530  printf("\n");
1531  }
1532  printf("\n");
1533  printf("K:\n");
1534  for (int i = 0; i < 6; ++i) {
1535  for (int j = 0; j < 2; ++j)
1536  printf("%8f ", K.At(0, i, j));
1537  printf("\n");
1538  }
1539  printf("\n");
1540  printf("outPar:\n");
1541  for (int i = 0; i < 6; ++i) {
1542  printf("%8f ", outPar.At(0, i, 0));
1543  }
1544  printf("\n");
1545  printf("outErr:\n");
1546  for (int i = 0; i < 6; ++i) {
1547  for (int j = 0; j < 6; ++j)
1548  printf("%8f ", outErr.At(0, i, j));
1549  printf("\n");
1550  }
1551  printf("\n");
1552  }
1553 #endif
1554  }
1555  }
Matriplex::Matriplex< float, LL, 2, NN > MPlexL2
Definition: Matrix.h:67
Definition: Electron.h:6
void squashPhiMPlex(MPlexLV &par, const int N_proc)
void invertCramerSym(MPlexSym< T, D, N > &A, double *determ=nullptr)
Definition: MatriplexSym.h:418
Matriplex::MatriplexSym< float, 2, NN > MPlex2S
Definition: Matrix.h:62
Matriplex::Matriplex< float, 2, 1, NN > MPlex2V
Definition: Matrix.h:61

◆ kalmanOperationPlane()

void mkfit::kalmanOperationPlane ( const int  kfOp,
const MPlexLS psErr,
const MPlexLV psPar,
const MPlexHS msErr,
const MPlexHV msPar,
const MPlexHV plNrm,
const MPlexHV plDir,
MPlexLS outErr,
MPlexLV outPar,
MPlexQF outChi2,
const int  N_proc 
)

Definition at line 1107 of file KalmanUtilsMPlex.cc.

References mps_fire::i, Matriplex::invertCramerSym(), dqmiolumiharvest::j, KFO_Calculate_Chi2, KFO_Local_Cov, KFO_Update_Params, dqmiodumpmetadata::n, NN, and squashPhiMPlex().

Referenced by kalmanComputeChi2Plane(), kalmanPropagateAndComputeChi2Plane(), kalmanPropagateAndUpdatePlane(), and kalmanUpdatePlane().

1117  {
1118 #ifdef DEBUG
1119  {
1120  dmutex_guard;
1121  printf("psPar:\n");
1122  for (int i = 0; i < 6; ++i) {
1123  printf("%8f ", psPar.constAt(0, 0, i));
1124  printf("\n");
1125  }
1126  printf("\n");
1127  printf("psErr:\n");
1128  for (int i = 0; i < 6; ++i) {
1129  for (int j = 0; j < 6; ++j)
1130  printf("%8f ", psErr.constAt(0, i, j));
1131  printf("\n");
1132  }
1133  printf("\n");
1134  printf("msPar:\n");
1135  for (int i = 0; i < 3; ++i) {
1136  printf("%8f ", msPar.constAt(0, 0, i));
1137  printf("\n");
1138  }
1139  printf("\n");
1140  printf("msErr:\n");
1141  for (int i = 0; i < 3; ++i) {
1142  for (int j = 0; j < 3; ++j)
1143  printf("%8f ", msErr.constAt(0, i, j));
1144  printf("\n");
1145  }
1146  printf("\n");
1147  }
1148 #endif
1149 
1150  // Rotate global point on tangent plane to cylinder
1151  // Tangent point is half way between hit and propagate position
1152 
1153  // Rotation matrix
1154  // D0 D1 D2
1155  // X0 X1 X2
1156  // N0 N1 N2
1157  // where D is the strip direction vector plDir, N is the normal plNrm, and X is the cross product between the two
1158 
1159  MPlex2H prj;
1160  for (int n = 0; n < NN; ++n) {
1161  prj(n, 0, 0) = plDir(n, 0, 0);
1162  prj(n, 0, 1) = plDir(n, 1, 0);
1163  prj(n, 0, 2) = plDir(n, 2, 0);
1164  prj(n, 1, 0) = plNrm(n, 1, 0) * plDir(n, 2, 0) - plNrm(n, 2, 0) * plDir(n, 1, 0);
1165  prj(n, 1, 1) = plNrm(n, 2, 0) * plDir(n, 0, 0) - plNrm(n, 0, 0) * plDir(n, 2, 0);
1166  prj(n, 1, 2) = plNrm(n, 0, 0) * plDir(n, 1, 0) - plNrm(n, 1, 0) * plDir(n, 0, 0);
1167  }
1168 
1169  MPlexHV res_glo; //position residual in global coordinates
1170  SubtractFirst3(msPar, psPar, res_glo);
1171 
1172  MPlexHS resErr_glo; //covariance sum in global position coordinates
1173  AddIntoUpperLeft3x3(psErr, msErr, resErr_glo);
1174 
1175  MPlex2V res_loc; //position residual in local coordinates
1176  RotateResidualsOnPlane(prj, res_glo, res_loc);
1177  MPlex2S resErr_loc; //covariance sum in local position coordinates
1178  MPlex2H temp2H;
1179  ProjectResErr(prj, resErr_glo, temp2H);
1180  ProjectResErrTransp(prj, temp2H, resErr_loc);
1181 
1182 #ifdef DEBUG
1183  {
1184  dmutex_guard;
1185  printf("prj:\n");
1186  for (int i = 0; i < 2; ++i) {
1187  for (int j = 0; j < 3; ++j)
1188  printf("%8f ", prj.At(0, i, j));
1189  printf("\n");
1190  }
1191  printf("\n");
1192  printf("res_glo:\n");
1193  for (int i = 0; i < 3; ++i) {
1194  printf("%8f ", res_glo.At(0, i, 0));
1195  }
1196  printf("\n");
1197  printf("resErr_glo:\n");
1198  for (int i = 0; i < 3; ++i) {
1199  for (int j = 0; j < 3; ++j)
1200  printf("%8f ", resErr_glo.At(0, i, j));
1201  printf("\n");
1202  }
1203  printf("\n");
1204  printf("res_loc:\n");
1205  for (int i = 0; i < 2; ++i) {
1206  printf("%8f ", res_loc.At(0, i, 0));
1207  }
1208  printf("\n");
1209  printf("temp2H:\n");
1210  for (int i = 0; i < 2; ++i) {
1211  for (int j = 0; j < 3; ++j)
1212  printf("%8f ", temp2H.At(0, i, j));
1213  printf("\n");
1214  }
1215  printf("\n");
1216  printf("resErr_loc:\n");
1217  for (int i = 0; i < 2; ++i) {
1218  for (int j = 0; j < 2; ++j)
1219  printf("%8f ", resErr_loc.At(0, i, j));
1220  printf("\n");
1221  }
1222  printf("\n");
1223  }
1224 #endif
1225 
1226  //invert the 2x2 matrix
1227  Matriplex::invertCramerSym(resErr_loc);
1228 
1229  if (kfOp & KFO_Calculate_Chi2) {
1230  Chi2Similarity(res_loc, resErr_loc, outChi2);
1231 
1232 #ifdef DEBUG
1233  {
1234  dmutex_guard;
1235  printf("resErr_loc (Inv):\n");
1236  for (int i = 0; i < 2; ++i) {
1237  for (int j = 0; j < 2; ++j)
1238  printf("%8f ", resErr_loc.At(0, i, j));
1239  printf("\n");
1240  }
1241  printf("\n");
1242  printf("chi2: %8f\n", outChi2.At(0, 0, 0));
1243  }
1244 #endif
1245  }
1246 
1247  if (kfOp & KFO_Update_Params) {
1248  MPlexLS psErrLoc = psErr;
1249 
1250  MPlexH2 tempH2;
1251  MPlexL2 K; // kalman gain, fixme should be L2
1252  KalmanHTG(prj, resErr_loc, tempH2); // intermediate term to get kalman gain (H^T*G)
1253  KalmanGain(psErrLoc, tempH2, K);
1254 
1255  MultResidualsAdd(K, psPar, res_loc, outPar);
1256 
1257  squashPhiMPlex(outPar, N_proc); // ensure phi is between |pi|
1258 
1259  MPlexLL tempLL;
1260  KHMult(K, prj, tempLL);
1261  KHC(tempLL, psErrLoc, outErr);
1262  outErr.subtract(psErrLoc, outErr);
1263 
1264 #ifdef DEBUG
1265  {
1266  dmutex_guard;
1267  if (kfOp & KFO_Local_Cov) {
1268  printf("psErrLoc:\n");
1269  for (int i = 0; i < 6; ++i) {
1270  for (int j = 0; j < 6; ++j)
1271  printf("% 8e ", psErrLoc.At(0, i, j));
1272  printf("\n");
1273  }
1274  printf("\n");
1275  }
1276  printf("resErr_loc (Inv):\n");
1277  for (int i = 0; i < 2; ++i) {
1278  for (int j = 0; j < 2; ++j)
1279  printf("%8f ", resErr_loc.At(0, i, j));
1280  printf("\n");
1281  }
1282  printf("\n");
1283  printf("tempH2:\n");
1284  for (int i = 0; i < 3; ++i) {
1285  for (int j = 0; j < 2; ++j)
1286  printf("%8f ", tempH2.At(0, i, j));
1287  printf("\n");
1288  }
1289  printf("\n");
1290  printf("K:\n");
1291  for (int i = 0; i < 6; ++i) {
1292  for (int j = 0; j < 2; ++j)
1293  printf("%8f ", K.At(0, i, j));
1294  printf("\n");
1295  }
1296  printf("\n");
1297  printf("tempLL:\n");
1298  for (int i = 0; i < 6; ++i) {
1299  for (int j = 0; j < 6; ++j)
1300  printf("%8f ", tempLL.At(0, i, j));
1301  printf("\n");
1302  }
1303  printf("\n");
1304  printf("outPar:\n");
1305  for (int i = 0; i < 6; ++i) {
1306  printf("%8f ", outPar.At(0, i, 0));
1307  }
1308  printf("\n");
1309  printf("outErr:\n");
1310  for (int i = 0; i < 6; ++i) {
1311  for (int j = 0; j < 6; ++j)
1312  printf("%8f ", outErr.At(0, i, j));
1313  printf("\n");
1314  }
1315  printf("\n");
1316  }
1317 #endif
1318  }
1319  }
Matriplex::Matriplex< float, LL, LL, NN > MPlexLL
Definition: Matrix.h:48
Matriplex::Matriplex< float, HH, 1, NN > MPlexHV
Definition: Matrix.h:53
Matriplex::Matriplex< float, LL, 2, NN > MPlexL2
Definition: Matrix.h:67
void squashPhiMPlex(MPlexLV &par, const int N_proc)
void invertCramerSym(MPlexSym< T, D, N > &A, double *determ=nullptr)
Definition: MatriplexSym.h:418
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
Matriplex::Matriplex< float, HH, 2, NN > MPlexH2
Definition: Matrix.h:68
Matriplex::Matriplex< float, 2, HH, NN > MPlex2H
Definition: Matrix.h:69
Matriplex::MatriplexSym< float, 2, NN > MPlex2S
Definition: Matrix.h:62
Matriplex::MatriplexSym< float, HH, NN > MPlexHS
Definition: Matrix.h:54
Matriplex::Matriplex< float, 2, 1, NN > MPlex2V
Definition: Matrix.h:61
Matriplex::MatriplexSym< float, LL, NN > MPlexLS
Definition: Matrix.h:50

◆ kalmanPropagateAndComputeChi2()

void mkfit::kalmanPropagateAndComputeChi2 ( const MPlexLS psErr,
const MPlexLV psPar,
const MPlexQI inChg,
const MPlexHS msErr,
const MPlexHV msPar,
MPlexQF outChi2,
MPlexLV propPar,
MPlexQI outFailFlag,
const int  N_proc,
const PropagationFlags propFlags,
const bool  propToHit 
)

Definition at line 763 of file KalmanUtilsMPlex.cc.

References Matriplex::hypot(), kalmanOperation(), KFO_Calculate_Chi2, dqmiodumpmetadata::n, NN, and propagateHelixToRMPlex().

773  {
774  propPar = psPar;
775  if (propToHit) {
776  MPlexLS propErr;
777  MPlexQF msRad;
778 #pragma omp simd
779  for (int n = 0; n < NN; ++n) {
780  msRad.At(n, 0, 0) = std::hypot(msPar.constAt(n, 0, 0), msPar.constAt(n, 1, 0));
781  }
782 
783  propagateHelixToRMPlex(psErr, psPar, inChg, msRad, propErr, propPar, outFailFlag, N_proc, propFlags);
784 
785  kalmanOperation(KFO_Calculate_Chi2, propErr, propPar, msErr, msPar, dummy_err, dummy_par, outChi2, N_proc);
786  } else {
787  kalmanOperation(KFO_Calculate_Chi2, psErr, psPar, msErr, msPar, dummy_err, dummy_par, outChi2, N_proc);
788  }
789  }
MPlex< T, D1, D2, N > hypot(const MPlex< T, D1, D2, N > &a, const MPlex< T, D1, D2, N > &b)
Definition: Matriplex.h:436
void propagateHelixToRMPlex(const MPlexLS &inErr, const MPlexLV &inPar, const MPlexQI &inChg, const MPlexQF &msRad, MPlexLS &outErr, MPlexLV &outPar, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &pflags, const MPlexQI *noMatEffPtr)
void kalmanOperation(const int kfOp, const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, MPlexLS &outErr, MPlexLV &outPar, MPlexQF &outChi2, const int N_proc)
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
Matriplex::Matriplex< float, 1, 1, NN > MPlexQF
Definition: Matrix.h:71
Matriplex::MatriplexSym< float, LL, NN > MPlexLS
Definition: Matrix.h:50

◆ kalmanPropagateAndComputeChi2Endcap()

void mkfit::kalmanPropagateAndComputeChi2Endcap ( const MPlexLS psErr,
const MPlexLV psPar,
const MPlexQI inChg,
const MPlexHS msErr,
const MPlexHV msPar,
MPlexQF outChi2,
MPlexLV propPar,
MPlexQI outFailFlag,
const int  N_proc,
const PropagationFlags propFlags,
const bool  propToHit 
)

Definition at line 1381 of file KalmanUtilsMPlex.cc.

References kalmanOperationEndcap(), KFO_Calculate_Chi2, dqmiodumpmetadata::n, NN, and propagateHelixToZMPlex().

1391  {
1392  propPar = psPar;
1393  if (propToHit) {
1394  MPlexLS propErr;
1395  MPlexQF msZ;
1396 #pragma omp simd
1397  for (int n = 0; n < NN; ++n) {
1398  msZ.At(n, 0, 0) = msPar.constAt(n, 2, 0);
1399  }
1400 
1401  propagateHelixToZMPlex(psErr, psPar, inChg, msZ, propErr, propPar, outFailFlag, N_proc, propFlags);
1402 
1403  kalmanOperationEndcap(KFO_Calculate_Chi2, propErr, propPar, msErr, msPar, dummy_err, dummy_par, outChi2, N_proc);
1404  } else {
1405  kalmanOperationEndcap(KFO_Calculate_Chi2, psErr, psPar, msErr, msPar, dummy_err, dummy_par, outChi2, N_proc);
1406  }
1407  }
void kalmanOperationEndcap(const int kfOp, const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, MPlexLS &outErr, MPlexLV &outPar, MPlexQF &outChi2, const int N_proc)
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
Matriplex::Matriplex< float, 1, 1, NN > MPlexQF
Definition: Matrix.h:71
void propagateHelixToZMPlex(const MPlexLS &inErr, const MPlexLV &inPar, const MPlexQI &inChg, const MPlexQF &msZ, MPlexLS &outErr, MPlexLV &outPar, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &pflags, const MPlexQI *noMatEffPtr)
Matriplex::MatriplexSym< float, LL, NN > MPlexLS
Definition: Matrix.h:50

◆ kalmanPropagateAndComputeChi2Plane()

void mkfit::kalmanPropagateAndComputeChi2Plane ( const MPlexLS psErr,
const MPlexLV psPar,
const MPlexQI inChg,
const MPlexHS msErr,
const MPlexHV msPar,
const MPlexHV plNrm,
const MPlexHV plDir,
MPlexQF outChi2,
MPlexLV propPar,
MPlexQI outFailFlag,
const int  N_proc,
const PropagationFlags propFlags,
const bool  propToHit 
)

Definition at line 1079 of file KalmanUtilsMPlex.cc.

References kalmanOperationPlane(), KFO_Calculate_Chi2, and propagateHelixToPlaneMPlex().

Referenced by mkfit::MkFinder::findCandidates(), and mkfit::MkFinder::findCandidatesCloneEngine().

1091  {
1092  propPar = psPar;
1093  if (propToHit) {
1094  MPlexLS propErr;
1095  propagateHelixToPlaneMPlex(psErr, psPar, inChg, msPar, plNrm, propErr, propPar, outFailFlag, N_proc, propFlags);
1096 
1098  KFO_Calculate_Chi2, propErr, propPar, msErr, msPar, plNrm, plDir, dummy_err, dummy_par, outChi2, N_proc);
1099  } else {
1101  KFO_Calculate_Chi2, psErr, psPar, msErr, msPar, plNrm, plDir, dummy_err, dummy_par, outChi2, N_proc);
1102  }
1103  }
void propagateHelixToPlaneMPlex(const MPlexLS &inErr, const MPlexLV &inPar, const MPlexQI &inChg, const MPlexHV &plPnt, const MPlexHV &plNrm, MPlexLS &outErr, MPlexLV &outPar, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &pflags, const MPlexQI *noMatEffPtr)
void kalmanOperationPlane(const int kfOp, const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, const MPlexHV &plNrm, const MPlexHV &plDir, MPlexLS &outErr, MPlexLV &outPar, MPlexQF &outChi2, const int N_proc)
Matriplex::MatriplexSym< float, LL, NN > MPlexLS
Definition: Matrix.h:50

◆ kalmanPropagateAndUpdate()

void mkfit::kalmanPropagateAndUpdate ( const MPlexLS psErr,
const MPlexLV psPar,
MPlexQI Chg,
const MPlexHS msErr,
const MPlexHV msPar,
MPlexLS outErr,
MPlexLV outPar,
MPlexQI outFailFlag,
const int  N_proc,
const PropagationFlags propFlags,
const bool  propToHit 
)

Definition at line 715 of file KalmanUtilsMPlex.cc.

References Matriplex::hypot(), kalmanOperation(), KFO_Local_Cov, KFO_Update_Params, dqmiodumpmetadata::n, NN, and propagateHelixToRMPlex().

725  {
726  if (propToHit) {
727  MPlexLS propErr;
728  MPlexLV propPar;
729  MPlexQF msRad;
730 #pragma omp simd
731  for (int n = 0; n < NN; ++n) {
732  msRad.At(n, 0, 0) = std::hypot(msPar.constAt(n, 0, 0), msPar.constAt(n, 1, 0));
733  }
734 
735  propagateHelixToRMPlex(psErr, psPar, Chg, msRad, propErr, propPar, outFailFlag, N_proc, propFlags);
736 
738  KFO_Update_Params | KFO_Local_Cov, propErr, propPar, msErr, msPar, outErr, outPar, dummy_chi2, N_proc);
739  } else {
741  KFO_Update_Params | KFO_Local_Cov, psErr, psPar, msErr, msPar, outErr, outPar, dummy_chi2, N_proc);
742  }
743  for (int n = 0; n < NN; ++n) {
744  if (outPar.At(n, 3, 0) < 0) {
745  Chg.At(n, 0, 0) = -Chg.At(n, 0, 0);
746  outPar.At(n, 3, 0) = -outPar.At(n, 3, 0);
747  }
748  }
749  }
MPlex< T, D1, D2, N > hypot(const MPlex< T, D1, D2, N > &a, const MPlex< T, D1, D2, N > &b)
Definition: Matriplex.h:436
void propagateHelixToRMPlex(const MPlexLS &inErr, const MPlexLV &inPar, const MPlexQI &inChg, const MPlexQF &msRad, MPlexLS &outErr, MPlexLV &outPar, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &pflags, const MPlexQI *noMatEffPtr)
Matriplex::Matriplex< float, LL, 1, NN > MPlexLV
Definition: Matrix.h:49
void kalmanOperation(const int kfOp, const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, MPlexLS &outErr, MPlexLV &outPar, MPlexQF &outChi2, const int N_proc)
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
Matriplex::Matriplex< float, 1, 1, NN > MPlexQF
Definition: Matrix.h:71
Matriplex::MatriplexSym< float, LL, NN > MPlexLS
Definition: Matrix.h:50

◆ kalmanPropagateAndUpdateEndcap()

void mkfit::kalmanPropagateAndUpdateEndcap ( const MPlexLS psErr,
const MPlexLV psPar,
MPlexQI Chg,
const MPlexHS msErr,
const MPlexHV msPar,
MPlexLS outErr,
MPlexLV outPar,
MPlexQI outFailFlag,
const int  N_proc,
const PropagationFlags propFlags,
const bool  propToHit 
)

Definition at line 1335 of file KalmanUtilsMPlex.cc.

References kalmanOperationEndcap(), KFO_Update_Params, dqmiodumpmetadata::n, NN, and propagateHelixToZMPlex().

1345  {
1346  if (propToHit) {
1347  MPlexLS propErr;
1348  MPlexLV propPar;
1349  MPlexQF msZ;
1350 #pragma omp simd
1351  for (int n = 0; n < NN; ++n) {
1352  msZ.At(n, 0, 0) = msPar.constAt(n, 2, 0);
1353  }
1354 
1355  propagateHelixToZMPlex(psErr, psPar, Chg, msZ, propErr, propPar, outFailFlag, N_proc, propFlags);
1356 
1357  kalmanOperationEndcap(KFO_Update_Params, propErr, propPar, msErr, msPar, outErr, outPar, dummy_chi2, N_proc);
1358  } else {
1359  kalmanOperationEndcap(KFO_Update_Params, psErr, psPar, msErr, msPar, outErr, outPar, dummy_chi2, N_proc);
1360  }
1361  for (int n = 0; n < NN; ++n) {
1362  if (outPar.At(n, 3, 0) < 0) {
1363  Chg.At(n, 0, 0) = -Chg.At(n, 0, 0);
1364  outPar.At(n, 3, 0) = -outPar.At(n, 3, 0);
1365  }
1366  }
1367  }
void kalmanOperationEndcap(const int kfOp, const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, MPlexLS &outErr, MPlexLV &outPar, MPlexQF &outChi2, const int N_proc)
Matriplex::Matriplex< float, LL, 1, NN > MPlexLV
Definition: Matrix.h:49
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
Matriplex::Matriplex< float, 1, 1, NN > MPlexQF
Definition: Matrix.h:71
void propagateHelixToZMPlex(const MPlexLS &inErr, const MPlexLV &inPar, const MPlexQI &inChg, const MPlexQF &msZ, MPlexLS &outErr, MPlexLV &outPar, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &pflags, const MPlexQI *noMatEffPtr)
Matriplex::MatriplexSym< float, LL, NN > MPlexLS
Definition: Matrix.h:50

◆ kalmanPropagateAndUpdatePlane()

void mkfit::kalmanPropagateAndUpdatePlane ( const MPlexLS psErr,
const MPlexLV psPar,
MPlexQI Chg,
const MPlexHS msErr,
const MPlexHV msPar,
const MPlexHV plNrm,
const MPlexHV plDir,
MPlexLS outErr,
MPlexLV outPar,
MPlexQI outFailFlag,
const int  N_proc,
const PropagationFlags propFlags,
const bool  propToHit 
)

Definition at line 1014 of file KalmanUtilsMPlex.cc.

References kalmanOperationPlane(), KFO_Local_Cov, KFO_Update_Params, dqmiodumpmetadata::n, NN, and propagateHelixToPlaneMPlex().

Referenced by mkfit::MkFinder::updateWithLoadedHit().

1026  {
1027  if (propToHit) {
1028  MPlexLS propErr;
1029  MPlexLV propPar;
1030  propagateHelixToPlaneMPlex(psErr, psPar, Chg, msPar, plNrm, propErr, propPar, outFailFlag, N_proc, propFlags);
1031 
1033  propErr,
1034  propPar,
1035  msErr,
1036  msPar,
1037  plNrm,
1038  plDir,
1039  outErr,
1040  outPar,
1041  dummy_chi2,
1042  N_proc);
1043  } else {
1045  psErr,
1046  psPar,
1047  msErr,
1048  msPar,
1049  plNrm,
1050  plDir,
1051  outErr,
1052  outPar,
1053  dummy_chi2,
1054  N_proc);
1055  }
1056  for (int n = 0; n < NN; ++n) {
1057  if (outPar.At(n, 3, 0) < 0) {
1058  Chg.At(n, 0, 0) = -Chg.At(n, 0, 0);
1059  outPar.At(n, 3, 0) = -outPar.At(n, 3, 0);
1060  }
1061  }
1062  }
void propagateHelixToPlaneMPlex(const MPlexLS &inErr, const MPlexLV &inPar, const MPlexQI &inChg, const MPlexHV &plPnt, const MPlexHV &plNrm, MPlexLS &outErr, MPlexLV &outPar, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &pflags, const MPlexQI *noMatEffPtr)
void kalmanOperationPlane(const int kfOp, const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, const MPlexHV &plNrm, const MPlexHV &plDir, MPlexLS &outErr, MPlexLV &outPar, MPlexQF &outChi2, const int N_proc)
Matriplex::Matriplex< float, LL, 1, NN > MPlexLV
Definition: Matrix.h:49
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
Matriplex::MatriplexSym< float, LL, NN > MPlexLS
Definition: Matrix.h:50

◆ kalmanUpdate()

void mkfit::kalmanUpdate ( const MPlexLS psErr,
const MPlexLV psPar,
const MPlexHS msErr,
const MPlexHV msPar,
MPlexLS outErr,
MPlexLV outPar,
const int  N_proc 
)

Definition at line 705 of file KalmanUtilsMPlex.cc.

References kalmanOperation(), KFO_Local_Cov, and KFO_Update_Params.

Referenced by mkfit::MkFitter::fitTracksWithInterSlurp().

711  {
712  kalmanOperation(KFO_Update_Params | KFO_Local_Cov, psErr, psPar, msErr, msPar, outErr, outPar, dummy_chi2, N_proc);
713  }
void kalmanOperation(const int kfOp, const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, MPlexLS &outErr, MPlexLV &outPar, MPlexQF &outChi2, const int N_proc)

◆ kalmanUpdateEndcap()

void mkfit::kalmanUpdateEndcap ( const MPlexLS psErr,
const MPlexLV psPar,
const MPlexHS msErr,
const MPlexHV msPar,
MPlexLS outErr,
MPlexLV outPar,
const int  N_proc 
)

Definition at line 1325 of file KalmanUtilsMPlex.cc.

References kalmanOperationEndcap(), and KFO_Update_Params.

1331  {
1332  kalmanOperationEndcap(KFO_Update_Params, psErr, psPar, msErr, msPar, outErr, outPar, dummy_chi2, N_proc);
1333  }
void kalmanOperationEndcap(const int kfOp, const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, MPlexLS &outErr, MPlexLV &outPar, MPlexQF &outChi2, const int N_proc)

◆ kalmanUpdatePlane()

void mkfit::kalmanUpdatePlane ( const MPlexLS psErr,
const MPlexLV psPar,
const MPlexHS msErr,
const MPlexHV msPar,
const MPlexHV plNrm,
const MPlexHV plDir,
MPlexLS outErr,
MPlexLV outPar,
const int  N_proc 
)

Definition at line 1001 of file KalmanUtilsMPlex.cc.

References kalmanOperationPlane(), KFO_Local_Cov, and KFO_Update_Params.

1009  {
1011  KFO_Update_Params | KFO_Local_Cov, psErr, psPar, msErr, msPar, plNrm, plDir, outErr, outPar, dummy_chi2, N_proc);
1012  }
void kalmanOperationPlane(const int kfOp, const MPlexLS &psErr, const MPlexLV &psPar, const MPlexHS &msErr, const MPlexHV &msPar, const MPlexHV &plNrm, const MPlexHV &plDir, MPlexLS &outErr, MPlexLV &outPar, MPlexQF &outChi2, const int N_proc)

◆ operator==()

template<class T , class U >
bool mkfit::operator== ( const CcAlloc< T > &  a,
const CcAlloc< U > &  b 
)

Definition at line 151 of file TrackStructures.h.

References a, and b.

151  {
152  return a.pool_id() == b.pool_id();
153  }
double b
Definition: hdecay.h:120
double a
Definition: hdecay.h:121

◆ passStripChargePCMfromTrack()

bool mkfit::passStripChargePCMfromTrack ( int  itrack,
bool  isBarrel,
unsigned int  pcm,
unsigned int  pcmMin,
const MPlexLV pPar,
const MPlexHS msErr 
)

Definition at line 1217 of file MkFinder.cc.

References funct::abs(), funct::cos(), dprint, f, PixelPluginsPhase0_cfi::isBarrel, mkfit::Hit::maxChargePerCM(), amptDefault_cfi::proj, funct::sin(), and mathSSE::sqrt().

Referenced by mkfit::MkFinder::findCandidates(), and mkfit::MkFinder::findCandidatesCloneEngine().

1218  {
1219  //skip the overflow case
1220  if (pcm >= Hit::maxChargePerCM())
1221  return true;
1222 
1223  float qSF;
1224  if (isBarrel) { //project in x,y, assuming zero-error direction is in this plane
1225  const float hitT2 = msErr.constAt(itrack, 0, 0) + msErr.constAt(itrack, 1, 1);
1226  const float hitT2inv = 1.f / hitT2;
1227  const float proj[3] = {msErr.constAt(itrack, 0, 0) * hitT2inv,
1228  msErr.constAt(itrack, 0, 1) * hitT2inv,
1229  msErr.constAt(itrack, 1, 1) * hitT2inv};
1230  const bool detXY_OK =
1231  std::abs(proj[0] * proj[2] - proj[1] * proj[1]) < 0.1f; //check that zero-direction is close
1232  const float cosP = cos(pPar.constAt(itrack, 4, 0));
1233  const float sinP = sin(pPar.constAt(itrack, 4, 0));
1234  const float sinT = std::abs(sin(pPar.constAt(itrack, 5, 0)));
1235  //qSF = sqrt[(px,py)*(1-proj)*(px,py)]/p = sinT*sqrt[(cosP,sinP)*(1-proj)*(cosP,sinP)].
1236  qSF = detXY_OK ? sinT * std::sqrt(std::abs(1.f + cosP * cosP * proj[0] + sinP * sinP * proj[2] -
1237  2.f * cosP * sinP * proj[1]))
1238  : 1.f;
1239  } else { //project on z
1240  // p_zLocal/p = p_z/p = cosT
1241  qSF = std::abs(cos(pPar.constAt(itrack, 5, 0)));
1242  }
1243 
1244  const float qCorr = pcm * qSF;
1245  dprint("pcm " << pcm << " * " << qSF << " = " << qCorr << " vs " << pcmMin);
1246  return qCorr > pcmMin;
1247  }
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
T sqrt(T t)
Definition: SSEVec.h:19
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
double f[11][100]
#define dprint(x)
Definition: Debug.h:95

◆ print() [1/5]

void mkfit::print ( std::string_view  label,
const MeasurementState s 
)

Definition at line 8 of file Hit.cc.

References gather_cfg::cout, dumpMatrix(), label, and alignCSCRings::s.

Referenced by mkfit::Shell::Compare(), conformalFitMPlex(), print(), and mkfit::StdSeq::track_print().

8  {
9  std::cout << label << std::endl;
10  std::cout << "x: " << s.parameters()[0] << " y: " << s.parameters()[1] << " z: " << s.parameters()[2] << std::endl
11  << "errors: " << std::endl;
12  dumpMatrix(s.errors());
13  std::cout << std::endl;
14  }
char const * label
void dumpMatrix(Matrix m)
Definition: MatrixSTypes.h:37

◆ print() [2/5]

void mkfit::print ( const TrackState s)

Definition at line 402 of file Track.cc.

References gather_cfg::cout, dumpMatrix(), and alignCSCRings::s.

402  {
403  std::cout << " x: " << s.parameters[0] << " y: " << s.parameters[1] << " z: " << s.parameters[2] << std::endl
404  << " px: " << s.parameters[3] << " py: " << s.parameters[4] << " pz: " << s.parameters[5] << std::endl
405  << "valid: " << s.valid << " errors: " << std::endl;
406  dumpMatrix(s.errors);
407  std::cout << std::endl;
408  }
void dumpMatrix(Matrix m)
Definition: MatrixSTypes.h:37

◆ print() [3/5]

void mkfit::print ( std::string  pfx,
int  itrack,
const Track trk,
bool  print_hits = false 
)

Definition at line 410 of file Track.cc.

References gather_cfg::cout, mkfit::Track::getHitIdx(), mkfit::Track::getHitLyr(), mps_fire::i, mkfit::TrackBase::label(), mkfit::Track::nFoundHits(), mkfit::Track::nTotalHits(), geometryDiffVisualizer::pfx, print(), and mkfit::TrackBase::state().

410  {
411  std::cout << std::endl
412  << pfx << ": " << itrack << " hits: " << trk.nFoundHits() << " label: " << trk.label() << " State"
413  << std::endl;
414  print(trk.state());
415  if (print_hits) {
416  for (int i = 0; i < trk.nTotalHits(); ++i)
417  printf(" %2d: lyr %2d idx %d\n", i, trk.getHitLyr(i), trk.getHitIdx(i));
418  }
419  }
void print(TMatrixD &m, const char *label=nullptr, bool mathematicaFormat=false)
Definition: Utilities.cc:47

◆ print() [4/5]

void mkfit::print ( std::string  pfx,
const TrackState s 
)

Definition at line 421 of file Track.cc.

References gather_cfg::cout, geometryDiffVisualizer::pfx, print(), and alignCSCRings::s.

421  {
422  std::cout << pfx << std::endl;
423  print(s);
424  }
void print(TMatrixD &m, const char *label=nullptr, bool mathematicaFormat=false)
Definition: Utilities.cc:47

◆ print() [5/5]

void mkfit::print ( std::string  pfx,
int  itrack,
const Track trk,
const Event ev 
)

Definition at line 1065 of file Event.cc.

References gather_cfg::cout, makeMEIFBenchmarkPlots::ev, mkfit::Track::getHitOnTrack(), h, mps_fire::i, mkfit::TrackBase::label(), mkfit::Track::nFoundHits(), mkfit::Track::nTotalHits(), geometryDiffVisualizer::pfx, print(), and mkfit::TrackBase::state().

1065  {
1066  std::cout << std::endl
1067  << pfx << ": " << itrack << " hits: " << trk.nFoundHits() << " label: " << trk.label()
1068  << " State:" << std::endl;
1069  print(trk.state());
1070 
1071  for (int i = 0; i < trk.nTotalHits(); ++i) {
1072  auto hot = trk.getHitOnTrack(i);
1073  printf(" %2d: lyr %2d idx %5d", i, hot.layer, hot.index);
1074  if (hot.index >= 0) {
1075  auto &h = ev.layerHits_[hot.layer][hot.index];
1076  int hl = ev.simHitsInfo_[h.mcHitID()].mcTrackID_;
1077  printf(" %4d %8.3f %8.3f %8.3f r=%.3f\n", hl, h.x(), h.y(), h.z(), h.r());
1078  } else {
1079  printf("\n");
1080  }
1081  }
1082  }
void print(std::string pfx, int itrack, const Track &trk, const Event &ev)
Definition: Event.cc:1065
The Signals That Services Can Subscribe To This is based on ActivityRegistry h
Helper function to determine trigger accepts.
Definition: Activities.doc:4

◆ propagateHelixToPlaneMPlex()

void mkfit::propagateHelixToPlaneMPlex ( const MPlexLS inErr,
const MPlexLV inPar,
const MPlexQI inChg,
const MPlexHV plPnt,
const MPlexHV plNrm,
MPlexLS outErr,
MPlexLV outPar,
MPlexQI outFailFlag,
const int  N_proc,
const PropagationFlags pflags,
const MPlexQI noMatEffPtr 
)

Definition at line 1082 of file PropagationMPlex.cc.

References funct::abs(), mkfit::PropagationFlags::apply_material, applyMaterialEffects(), funct::cos(), debug, dprint_np, dprintf, g_debug, helixAtPlane(), Matriplex::hypot(), mps_fire::i, dqmiolumiharvest::j, GetRecoTauVFromDQM_MC_cff::kk, mkfit::TrackerInfo::material_checked(), dqmiodumpmetadata::n, NN, mkfit::TrackerInfo::Material::radl, funct::sin(), mathSSE::sqrt(), squashPhiMPlex(), funct::tan(), groupFilesInBlocks::temp, and mkfit::PropagationFlags::tracker_info.

Referenced by kalmanPropagateAndComputeChi2Plane(), and kalmanPropagateAndUpdatePlane().

1092  {
1093  // debug = true;
1094 
1095  outErr = inErr;
1096  outPar = inPar;
1097 
1098  MPlexQF pathL;
1099  MPlexLL errorProp;
1100 
1101  helixAtPlane(inPar, inChg, plPnt, plNrm, pathL, outPar, errorProp, outFailFlag, N_proc, pflags);
1102 
1103  for (int n = 0; n < NN; ++n) {
1104  dprint_np(
1105  n,
1106  "propagation to plane end, dump parameters\n"
1107  //<< " D = " << s[n] << " alpha = " << s[n] * std::sin(inPar(n, 5, 0)) * inPar(n, 3, 0) * kinv[n] << " kinv = " << kinv[n] << std::endl
1108  << " pos = " << outPar(n, 0, 0) << " " << outPar(n, 1, 0) << " " << outPar(n, 2, 0) << "\t\t r="
1109  << std::sqrt(outPar(n, 0, 0) * outPar(n, 0, 0) + outPar(n, 1, 0) * outPar(n, 1, 0)) << std::endl
1110  << " mom = " << outPar(n, 3, 0) << " " << outPar(n, 4, 0) << " " << outPar(n, 5, 0) << std::endl
1111  << " cart= " << std::cos(outPar(n, 4, 0)) / outPar(n, 3, 0) << " "
1112  << std::sin(outPar(n, 4, 0)) / outPar(n, 3, 0) << " " << 1. / (outPar(n, 3, 0) * tan(outPar(n, 5, 0)))
1113  << "\t\tpT=" << 1. / std::abs(outPar(n, 3, 0)) << std::endl);
1114  }
1115 
1116 #ifdef DEBUG
1117  if (debug && g_debug) {
1118  for (int kk = 0; kk < N_proc; ++kk) {
1119  dprintf("inPar %d\n", kk);
1120  for (int i = 0; i < 6; ++i) {
1121  dprintf("%8f ", inPar.constAt(kk, i, 0));
1122  }
1123  dprintf("\n");
1124  dprintf("inErr %d\n", kk);
1125  for (int i = 0; i < 6; ++i) {
1126  for (int j = 0; j < 6; ++j)
1127  dprintf("%8f ", inErr.constAt(kk, i, j));
1128  dprintf("\n");
1129  }
1130  dprintf("\n");
1131 
1132  for (int kk = 0; kk < N_proc; ++kk) {
1133  dprintf("plNrm %d\n", kk);
1134  for (int j = 0; j < 3; ++j)
1135  dprintf("%8f ", plNrm.constAt(kk, 0, j));
1136  }
1137  dprintf("\n");
1138 
1139  for (int kk = 0; kk < N_proc; ++kk) {
1140  dprintf("pathL %d\n", kk);
1141  for (int j = 0; j < 1; ++j)
1142  dprintf("%8f ", pathL.constAt(kk, 0, j));
1143  }
1144  dprintf("\n");
1145 
1146  dprintf("errorProp %d\n", kk);
1147  for (int i = 0; i < 6; ++i) {
1148  for (int j = 0; j < 6; ++j)
1149  dprintf("%8f ", errorProp.At(kk, i, j));
1150  dprintf("\n");
1151  }
1152  dprintf("\n");
1153  }
1154  }
1155 #endif
1156 
1157  // Matriplex version of:
1158  // result.errors = ROOT::Math::Similarity(errorProp, outErr);
1159  MPlexLL temp;
1160  MultHelixPropFull(errorProp, outErr, temp);
1161  MultHelixPropTranspFull(errorProp, temp, outErr);
1162 
1163 #ifdef DEBUG
1164  if (debug && g_debug) {
1165  for (int kk = 0; kk < N_proc; ++kk) {
1166  dprintf("outErr %d\n", kk);
1167  for (int i = 0; i < 6; ++i) {
1168  for (int j = 0; j < 6; ++j)
1169  dprintf("%8f ", outErr.constAt(kk, i, j));
1170  dprintf("\n");
1171  }
1172  dprintf("\n");
1173  }
1174  }
1175 #endif
1176 
1177  if (pflags.apply_material) {
1178  MPlexQF hitsRl;
1179  MPlexQF hitsXi;
1180  MPlexQF propSign;
1181 
1182  const TrackerInfo& tinfo = *pflags.tracker_info;
1183 
1184 #pragma omp simd
1185  for (int n = 0; n < NN; ++n) {
1186  if (n >= N_proc || (noMatEffPtr && noMatEffPtr->constAt(n, 0, 0))) {
1187  hitsRl(n, 0, 0) = 0.f;
1188  hitsXi(n, 0, 0) = 0.f;
1189  } else {
1190  const float hypo = std::hypot(outPar(n, 0, 0), outPar(n, 1, 0));
1191  auto mat = tinfo.material_checked(std::abs(outPar(n, 2, 0)), hypo);
1192  hitsRl(n, 0, 0) = mat.radl;
1193  hitsXi(n, 0, 0) = mat.bbxi;
1194  }
1195  propSign(n, 0, 0) = (pathL(n, 0, 0) > 0.f ? 1.f : -1.f);
1196  }
1197  applyMaterialEffects(hitsRl, hitsXi, propSign, plNrm, outErr, outPar, N_proc);
1198 #ifdef DEBUG
1199  if (debug && g_debug) {
1200  for (int kk = 0; kk < N_proc; ++kk) {
1201  dprintf("propSign %d\n", kk);
1202  for (int i = 0; i < 1; ++i) {
1203  dprintf("%8f ", propSign.constAt(kk, i, 0));
1204  }
1205  dprintf("\n");
1206  dprintf("plNrm %d\n", kk);
1207  for (int i = 0; i < 3; ++i) {
1208  dprintf("%8f ", plNrm.constAt(kk, i, 0));
1209  }
1210  dprintf("\n");
1211  dprintf("outErr(after material) %d\n", kk);
1212  for (int i = 0; i < 6; ++i) {
1213  for (int j = 0; j < 6; ++j)
1214  dprintf("%8f ", outErr.constAt(kk, i, j));
1215  dprintf("\n");
1216  }
1217  dprintf("\n");
1218  }
1219  }
1220 #endif
1221  }
1222 
1223  squashPhiMPlex(outPar, N_proc); // ensure phi is between |pi|
1224 
1225  // PROP-FAIL-ENABLE To keep physics changes minimal, we always restore the
1226  // state to input when propagation fails -- as was the default before.
1227  // if (pflags.copy_input_state_on_fail) {
1228  for (int i = 0; i < N_proc; ++i) {
1229  if (outFailFlag(i, 0, 0)) {
1230  outPar.copySlot(i, inPar);
1231  outErr.copySlot(i, inErr);
1232  }
1233  }
1234  // }
1235  }
Matriplex::Matriplex< float, LL, LL, NN > MPlexLL
Definition: Matrix.h:48
MPlex< T, D1, D2, N > hypot(const MPlex< T, D1, D2, N > &a, const MPlex< T, D1, D2, N > &b)
Definition: Matriplex.h:436
#define dprint_np(n, x)
Definition: Debug.h:96
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
const TrackerInfo * tracker_info
void applyMaterialEffects(const MPlexQF &hitsRl, const MPlexQF &hitsXi, const MPlexQF &propSign, const MPlexHV &plNrm, MPlexLS &outErr, MPlexLV &outPar, const int N_proc)
void squashPhiMPlex(MPlexLV &par, const int N_proc)
bool g_debug
Definition: Debug.cc:2
T sqrt(T t)
Definition: SSEVec.h:19
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
Tan< T >::type tan(const T &t)
Definition: Tan.h:22
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
#define debug
Definition: HDRShower.cc:19
Matriplex::Matriplex< float, 1, 1, NN > MPlexQF
Definition: Matrix.h:71
void helixAtPlane(const MPlexLV &inPar, const MPlexQI &inChg, const MPlexHV &plPnt, const MPlexHV &plNrm, MPlexQF &pathL, MPlexLV &outPar, MPlexLL &errorProp, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &pflags)
Material material_checked(float z, float r) const
Definition: TrackerInfo.h:237
#define dprintf(...)
Definition: Debug.h:98

◆ propagateHelixToRMPlex()

void mkfit::propagateHelixToRMPlex ( const MPlexLS inErr,
const MPlexLV inPar,
const MPlexQI inChg,
const MPlexQF msRad,
MPlexLS outErr,
MPlexLV outPar,
MPlexQI outFailFlag,
const int  N_proc,
const PropagationFlags pflags,
const MPlexQI noMatEffPtr 
)

Definition at line 494 of file PropagationMPlex.cc.

References funct::abs(), mkfit::PropagationFlags::apply_material, applyMaterialEffects(), funct::cos(), debug, dprintf, g_debug, helixAtRFromIterativeCCS(), hipo(), mps_fire::i, dqmiolumiharvest::j, GetRecoTauVFromDQM_MC_cff::kk, mkfit::TrackerInfo::material_checked(), dqmiodumpmetadata::n, NN, mkfit::TrackerInfo::Material::radl, funct::sin(), squashPhiMPlex(), groupFilesInBlocks::temp, and mkfit::PropagationFlags::tracker_info.

Referenced by kalmanPropagateAndComputeChi2(), kalmanPropagateAndUpdate(), mkfit::MkBase::propagateTracksToHitR(), and mkfit::MkBase::propagateTracksToR().

503  {
504  // bool debug = true;
505 
506  // This is used further down when calculating similarity with errorProp (and before in DEBUG).
507  // MT: I don't think this really needed if we use inErr where required.
508  outErr = inErr;
509  // This requirement for helixAtRFromIterativeCCS_impl() and for helixAtRFromIterativeCCSFullJac().
510  // MT: This should be properly handled in both functions (expecting input in output parameters sucks).
511  outPar = inPar;
512 
513  MPlexLL errorProp;
514 
515  helixAtRFromIterativeCCS(inPar, inChg, msRad, outPar, errorProp, outFailFlag, N_proc, pflags);
516 
517 #ifdef DEBUG
518  if (debug && g_debug) {
519  for (int kk = 0; kk < N_proc; ++kk) {
520  dprintf("outErr before prop %d\n", kk);
521  for (int i = 0; i < 6; ++i) {
522  for (int j = 0; j < 6; ++j)
523  dprintf("%8f ", outErr.At(kk, i, j));
524  dprintf("\n");
525  }
526  dprintf("\n");
527 
528  dprintf("errorProp %d\n", kk);
529  for (int i = 0; i < 6; ++i) {
530  for (int j = 0; j < 6; ++j)
531  dprintf("%8f ", errorProp.At(kk, i, j));
532  dprintf("\n");
533  }
534  dprintf("\n");
535  }
536  }
537 #endif
538 
539  // MultHelixProp can be optimized for CCS coordinates, see GenMPlexOps.pl
540  MPlexLL temp;
541  MultHelixProp(errorProp, outErr, temp);
542  MultHelixPropTransp(errorProp, temp, outErr);
543  // can replace with: MultHelixPropFull(errorProp, outErr, temp); MultHelixPropTranspFull(errorProp, temp, outErr);
544 
545  if (pflags.apply_material) {
546  MPlexQF hitsRl;
547  MPlexQF hitsXi;
548  MPlexQF propSign;
549 
550  const TrackerInfo& tinfo = *pflags.tracker_info;
551 
552 #pragma omp simd
553  for (int n = 0; n < NN; ++n) {
554  if (n >= N_proc || (outFailFlag(n, 0, 0) || (noMatEffPtr && noMatEffPtr->constAt(n, 0, 0)))) {
555  hitsRl(n, 0, 0) = 0.f;
556  hitsXi(n, 0, 0) = 0.f;
557  } else {
558  auto mat = tinfo.material_checked(std::abs(outPar(n, 2, 0)), msRad(n, 0, 0));
559  hitsRl(n, 0, 0) = mat.radl;
560  hitsXi(n, 0, 0) = mat.bbxi;
561  }
562  const float r0 = hipo(inPar(n, 0, 0), inPar(n, 1, 0));
563  const float r = msRad(n, 0, 0);
564  propSign(n, 0, 0) = (r > r0 ? 1. : -1.);
565  }
566  MPlexHV plNrm;
567 #pragma omp simd
568  for (int n = 0; n < NN; ++n) {
569  plNrm(n, 0, 0) = std::cos(outPar.constAt(n, 4, 0));
570  plNrm(n, 1, 0) = std::sin(outPar.constAt(n, 4, 0));
571  plNrm(n, 2, 0) = 0.f;
572  }
573  applyMaterialEffects(hitsRl, hitsXi, propSign, plNrm, outErr, outPar, N_proc);
574  }
575 
576  squashPhiMPlex(outPar, N_proc); // ensure phi is between |pi|
577 
578  // Matriplex version of:
579  // result.errors = ROOT::Math::Similarity(errorProp, outErr);
580 
581  /*
582  // To be used with: MPT_DIM = 1
583  if (fabs(sqrt(outPar[0]*outPar[0]+outPar[1]*outPar[1]) - msRad[0]) > 0.0001)
584  {
585  std::cout << "DID NOT GET TO R, FailFlag=" << failFlag[0]
586  << " dR=" << msRad[0] - std::hypot(outPar[0],outPar[1])
587  << " r=" << msRad[0] << " rin=" << std::hypot(inPar[0],inPar[1]) << " rout=" << std::hypot(outPar[0],outPar[1])
588  << std::endl;
589  // std::cout << " pt=" << pt << " pz=" << inPar.At(n, 2) << std::endl;
590  }
591  */
592 
593  // PROP-FAIL-ENABLE To keep physics changes minimal, we always restore the
594  // state to input when propagation fails -- as was the default before.
595  // if (pflags.copy_input_state_on_fail) {
596  for (int i = 0; i < N_proc; ++i) {
597  if (outFailFlag(i, 0, 0)) {
598  outPar.copySlot(i, inPar);
599  outErr.copySlot(i, inErr);
600  }
601  }
602  // }
603  }
Matriplex::Matriplex< float, LL, LL, NN > MPlexLL
Definition: Matrix.h:48
Matriplex::Matriplex< float, HH, 1, NN > MPlexHV
Definition: Matrix.h:53
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
const TrackerInfo * tracker_info
void applyMaterialEffects(const MPlexQF &hitsRl, const MPlexQF &hitsXi, const MPlexQF &propSign, const MPlexHV &plNrm, MPlexLS &outErr, MPlexLV &outPar, const int N_proc)
void squashPhiMPlex(MPlexLV &par, const int N_proc)
bool g_debug
Definition: Debug.cc:2
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
#define debug
Definition: HDRShower.cc:19
Matriplex::Matriplex< float, 1, 1, NN > MPlexQF
Definition: Matrix.h:71
float hipo(float x, float y)
Definition: Matrix.h:9
Material material_checked(float z, float r) const
Definition: TrackerInfo.h:237
void helixAtRFromIterativeCCS(const MPlexLV &inPar, const MPlexQI &inChg, const MPlexQF &msRad, MPlexLV &outPar, MPlexLL &errorProp, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &pflags)
#define dprintf(...)
Definition: Debug.h:98

◆ propagateHelixToZMPlex()

void mkfit::propagateHelixToZMPlex ( const MPlexLS inErr,
const MPlexLV inPar,
const MPlexQI inChg,
const MPlexQF msZ,
MPlexLS outErr,
MPlexLV outPar,
MPlexQI outFailFlag,
const int  N_proc,
const PropagationFlags pflags,
const MPlexQI noMatEffPtr 
)

Definition at line 607 of file PropagationMPlex.cc.

References funct::abs(), mkfit::PropagationFlags::apply_material, applyMaterialEffects(), debug, dprintf, g_debug, helixAtZ(), Matriplex::hypot(), mps_fire::i, dqmiolumiharvest::j, GetRecoTauVFromDQM_MC_cff::kk, mkfit::TrackerInfo::material_checked(), dqmiodumpmetadata::n, NN, mkfit::TrackerInfo::Material::radl, squashPhiMPlex(), groupFilesInBlocks::temp, and mkfit::PropagationFlags::tracker_info.

Referenced by kalmanPropagateAndComputeChi2Endcap(), kalmanPropagateAndUpdateEndcap(), mkfit::MkBase::propagateTracksToHitZ(), mkfit::MkBase::propagateTracksToPCAZ(), and mkfit::MkBase::propagateTracksToZ().

616  {
617  // debug = true;
618 
619  outErr = inErr;
620  outPar = inPar;
621 
622  MPlexLL errorProp;
623 
624  //helixAtZ_new(inPar, inChg, msZ, outPar, errorProp, outFailFlag, N_proc, pflags);
625  helixAtZ(inPar, inChg, msZ, outPar, errorProp, outFailFlag, N_proc, pflags);
626 
627 #ifdef DEBUG
628  if (debug && g_debug) {
629  for (int kk = 0; kk < N_proc; ++kk) {
630  dprintf("inPar %d\n", kk);
631  for (int i = 0; i < 6; ++i) {
632  dprintf("%8f ", inPar.constAt(kk, i, 0));
633  }
634  dprintf("\n");
635 
636  dprintf("inErr %d\n", kk);
637  for (int i = 0; i < 6; ++i) {
638  for (int j = 0; j < 6; ++j)
639  dprintf("%8f ", inErr.constAt(kk, i, j));
640  dprintf("\n");
641  }
642  dprintf("\n");
643 
644  dprintf("errorProp %d\n", kk);
645  for (int i = 0; i < 6; ++i) {
646  for (int j = 0; j < 6; ++j)
647  dprintf("%8f ", errorProp.At(kk, i, j));
648  dprintf("\n");
649  }
650  dprintf("\n");
651  }
652  }
653 #endif
654 
655 #ifdef DEBUG
656  if (debug && g_debug) {
657  for (int kk = 0; kk < N_proc; ++kk) {
658  dprintf("outErr %d\n", kk);
659  for (int i = 0; i < 6; ++i) {
660  for (int j = 0; j < 6; ++j)
661  dprintf("%8f ", outErr.constAt(kk, i, j));
662  dprintf("\n");
663  }
664  dprintf("\n");
665  }
666  }
667 #endif
668 
669  // Matriplex version of: result.errors = ROOT::Math::Similarity(errorProp, outErr);
670  MPlexLL temp;
671  MultHelixPropEndcap(errorProp, outErr, temp);
672  MultHelixPropTranspEndcap(errorProp, temp, outErr);
673  // can replace with: MultHelixPropFull(errorProp, outErr, temp); MultHelixPropTranspFull(errorProp, temp, outErr);
674 
675  if (pflags.apply_material) {
676  MPlexQF hitsRl;
677  MPlexQF hitsXi;
678  MPlexQF propSign;
679 
680  const TrackerInfo& tinfo = *pflags.tracker_info;
681 
682 #pragma omp simd
683  for (int n = 0; n < NN; ++n) {
684  if (n >= N_proc || (noMatEffPtr && noMatEffPtr->constAt(n, 0, 0))) {
685  hitsRl(n, 0, 0) = 0.f;
686  hitsXi(n, 0, 0) = 0.f;
687  } else {
688  const float hypo = std::hypot(outPar(n, 0, 0), outPar(n, 1, 0));
689  auto mat = tinfo.material_checked(std::abs(msZ(n, 0, 0)), hypo);
690  hitsRl(n, 0, 0) = mat.radl;
691  hitsXi(n, 0, 0) = mat.bbxi;
692  }
693  const float zout = msZ.constAt(n, 0, 0);
694  const float zin = inPar.constAt(n, 2, 0);
695  propSign(n, 0, 0) = (std::abs(zout) > std::abs(zin) ? 1.f : -1.f);
696  }
697  MPlexHV plNrm;
698 #pragma omp simd
699  for (int n = 0; n < NN; ++n) {
700  plNrm(n, 0, 0) = 0.f;
701  plNrm(n, 1, 0) = 0.f;
702  plNrm(n, 2, 0) = 1.f;
703  }
704  applyMaterialEffects(hitsRl, hitsXi, propSign, plNrm, outErr, outPar, N_proc);
705 #ifdef DEBUG
706  if (debug && g_debug) {
707  for (int kk = 0; kk < N_proc; ++kk) {
708  dprintf("propSign %d\n", kk);
709  for (int i = 0; i < 1; ++i) {
710  dprintf("%8f ", propSign.constAt(kk, i, 0));
711  }
712  dprintf("\n");
713  dprintf("plNrm %d\n", kk);
714  for (int i = 0; i < 3; ++i) {
715  dprintf("%8f ", plNrm.constAt(kk, i, 0));
716  }
717  dprintf("\n");
718  dprintf("outErr(after material) %d\n", kk);
719  for (int i = 0; i < 6; ++i) {
720  for (int j = 0; j < 6; ++j)
721  dprintf("%8f ", outErr.constAt(kk, i, j));
722  dprintf("\n");
723  }
724  dprintf("\n");
725  }
726  }
727 #endif
728  }
729 
730  squashPhiMPlex(outPar, N_proc); // ensure phi is between |pi|
731 
732  // PROP-FAIL-ENABLE To keep physics changes minimal, we always restore the
733  // state to input when propagation fails -- as was the default before.
734  // if (pflags.copy_input_state_on_fail) {
735  for (int i = 0; i < N_proc; ++i) {
736  if (outFailFlag(i, 0, 0)) {
737  outPar.copySlot(i, inPar);
738  outErr.copySlot(i, inErr);
739  }
740  }
741  // }
742  }
Matriplex::Matriplex< float, LL, LL, NN > MPlexLL
Definition: Matrix.h:48
MPlex< T, D1, D2, N > hypot(const MPlex< T, D1, D2, N > &a, const MPlex< T, D1, D2, N > &b)
Definition: Matriplex.h:436
Matriplex::Matriplex< float, HH, 1, NN > MPlexHV
Definition: Matrix.h:53
const TrackerInfo * tracker_info
void applyMaterialEffects(const MPlexQF &hitsRl, const MPlexQF &hitsXi, const MPlexQF &propSign, const MPlexHV &plNrm, MPlexLS &outErr, MPlexLV &outPar, const int N_proc)
void squashPhiMPlex(MPlexLV &par, const int N_proc)
bool g_debug
Definition: Debug.cc:2
void helixAtZ(const MPlexLV &inPar, const MPlexQI &inChg, const MPlexQF &msZ, MPlexLV &outPar, MPlexLL &errorProp, MPlexQI &outFailFlag, const int N_proc, const PropagationFlags &pflags)
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
#define debug
Definition: HDRShower.cc:19
Matriplex::Matriplex< float, 1, 1, NN > MPlexQF
Definition: Matrix.h:71
Material material_checked(float z, float r) const
Definition: TrackerInfo.h:237
#define dprintf(...)
Definition: Debug.h:98

◆ propagateLineToRMPlex()

void mkfit::propagateLineToRMPlex ( const MPlexLS psErr,
const MPlexLV psPar,
const MPlexHS msErr,
const MPlexHV msPar,
MPlexLS outErr,
MPlexLV outPar,
const int  N_proc 
)

Definition at line 19 of file PropagationMPlex.cc.

References A, B, dprint_np, hipo(), N, dqmiodumpmetadata::n, NN, AlCaHLTBitMon_ParallelJobs::p, DiDispStaMuonMonitor_cfi::pt, and mathSSE::sqrt().

25  {
26  // XXX Regenerate parts below with a script.
27 
28  const Matriplex::idx_t N = NN;
29 
30 #pragma omp simd
31  for (int n = 0; n < NN; ++n) {
32  const float cosA = (psPar[0 * N + n] * psPar[3 * N + n] + psPar[1 * N + n] * psPar[4 * N + n]) /
33  (std::sqrt((psPar[0 * N + n] * psPar[0 * N + n] + psPar[1 * N + n] * psPar[1 * N + n]) *
34  (psPar[3 * N + n] * psPar[3 * N + n] + psPar[4 * N + n] * psPar[4 * N + n])));
35  const float dr = (hipo(msPar[0 * N + n], msPar[1 * N + n]) - hipo(psPar[0 * N + n], psPar[1 * N + n])) / cosA;
36 
37  dprint_np(n, "propagateLineToRMPlex dr=" << dr);
38 
39  const float pt = hipo(psPar[3 * N + n], psPar[4 * N + n]);
40  const float p = dr / pt; // path
41  const float psq = p * p;
42 
43  outPar[0 * N + n] = psPar[0 * N + n] + p * psPar[3 * N + n];
44  outPar[1 * N + n] = psPar[1 * N + n] + p * psPar[4 * N + n];
45  outPar[2 * N + n] = psPar[2 * N + n] + p * psPar[5 * N + n];
46  outPar[3 * N + n] = psPar[3 * N + n];
47  outPar[4 * N + n] = psPar[4 * N + n];
48  outPar[5 * N + n] = psPar[5 * N + n];
49 
50  {
51  const MPlexLS& A = psErr;
52  MPlexLS& B = outErr;
53 
54  B.fArray[0 * N + n] = A.fArray[0 * N + n];
55  B.fArray[1 * N + n] = A.fArray[1 * N + n];
56  B.fArray[2 * N + n] = A.fArray[2 * N + n];
57  B.fArray[3 * N + n] = A.fArray[3 * N + n];
58  B.fArray[4 * N + n] = A.fArray[4 * N + n];
59  B.fArray[5 * N + n] = A.fArray[5 * N + n];
60  B.fArray[6 * N + n] = A.fArray[6 * N + n] + p * A.fArray[0 * N + n];
61  B.fArray[7 * N + n] = A.fArray[7 * N + n] + p * A.fArray[1 * N + n];
62  B.fArray[8 * N + n] = A.fArray[8 * N + n] + p * A.fArray[3 * N + n];
63  B.fArray[9 * N + n] =
64  A.fArray[9 * N + n] + p * (A.fArray[6 * N + n] + A.fArray[6 * N + n]) + psq * A.fArray[0 * N + n];
65  B.fArray[10 * N + n] = A.fArray[10 * N + n] + p * A.fArray[1 * N + n];
66  B.fArray[11 * N + n] = A.fArray[11 * N + n] + p * A.fArray[2 * N + n];
67  B.fArray[12 * N + n] = A.fArray[12 * N + n] + p * A.fArray[4 * N + n];
68  B.fArray[13 * N + n] =
69  A.fArray[13 * N + n] + p * (A.fArray[7 * N + n] + A.fArray[10 * N + n]) + psq * A.fArray[1 * N + n];
70  B.fArray[14 * N + n] =
71  A.fArray[14 * N + n] + p * (A.fArray[11 * N + n] + A.fArray[11 * N + n]) + psq * A.fArray[2 * N + n];
72  B.fArray[15 * N + n] = A.fArray[15 * N + n] + p * A.fArray[3 * N + n];
73  B.fArray[16 * N + n] = A.fArray[16 * N + n] + p * A.fArray[4 * N + n];
74  B.fArray[17 * N + n] = A.fArray[17 * N + n] + p * A.fArray[5 * N + n];
75  B.fArray[18 * N + n] =
76  A.fArray[18 * N + n] + p * (A.fArray[8 * N + n] + A.fArray[15 * N + n]) + psq * A.fArray[3 * N + n];
77  B.fArray[19 * N + n] =
78  A.fArray[19 * N + n] + p * (A.fArray[12 * N + n] + A.fArray[16 * N + n]) + psq * A.fArray[4 * N + n];
79  B.fArray[20 * N + n] =
80  A.fArray[20 * N + n] + p * (A.fArray[17 * N + n] + A.fArray[17 * N + n]) + psq * A.fArray[5 * N + n];
81  }
82 
83  dprint_np(n, "propagateLineToRMPlex arrive at r=" << hipo(outPar[0 * N + n], outPar[1 * N + n]));
84  }
85  }
Definition: APVGainStruct.h:7
#define dprint_np(n, x)
Definition: Debug.h:96
T sqrt(T t)
Definition: SSEVec.h:19
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
#define N
Definition: blowfish.cc:9
float hipo(float x, float y)
Definition: Matrix.h:9
Matriplex::MatriplexSym< float, LL, NN > MPlexLS
Definition: Matrix.h:50
Definition: APVGainStruct.h:7

◆ run_OneIteration()

void mkfit::run_OneIteration ( const TrackerInfo trackerInfo,
const IterationConfig itconf,
const EventOfHits eoh,
const std::vector< const std::vector< bool > *> &  hit_masks,
MkBuilder builder,
TrackVec seeds,
TrackVec out_tracks,
bool  do_seed_clean,
bool  do_backward_fit,
bool  do_remove_duplicates 
)

Definition at line 29 of file runFunctions.cc.

References mkfit::MkBuilder::backwardFit(), mkfit::MkBuilder::begin_event(), mkfit::MkBuilder::beginBkwSearch(), mkfit::MkBuilder::compactifyHitStorageForBestCand(), mkfit::MkBuilder::end_event(), mkfit::MkBuilder::endBkwSearch(), mkfit::MkBuilder::export_best_comb_cands(), mkfit::MkBuilder::filter_comb_cands(), mkfit::MkBuilder::find_tracks_load_seeds(), mkfit::MkBuilder::findTracksCloneEngine(), mkfit::SteeringParams::IT_BkwSearch, mkfit::IterationConfig::m_backward_drop_seed_hits, mkfit::IterationConfig::m_backward_fit_min_hits, mkfit::IterationConfig::m_backward_search, mkfit::IterationConfig::m_duplicate_cleaner, mkfit::IterationConfig::m_post_bkfit_filter, mkfit::IterationConfig::m_pre_bkfit_filter, mkfit::IterationConfig::m_requires_seed_hit_sorting, mkfit::IterationConfig::m_seed_cleaner, mkfit::EventOfHits::refBeamSpot(), mkfit::MkBuilder::release_memory(), alignCSCRings::s, mkfit::MkBuilder::seed_post_cleaning(), and DetachedQuadStep_cff::seeds.

Referenced by MkFitProducer::produce().

38  {
39  IterationMaskIfcCmssw it_mask_ifc(trackerInfo, hit_masks);
40 
41  MkJob job({trackerInfo, itconf, eoh, eoh.refBeamSpot(), &it_mask_ifc});
42 
43  builder.begin_event(&job, nullptr, __func__);
44 
45  // Seed cleaning not done on all iterations.
46  do_seed_clean = do_seed_clean && itconf.m_seed_cleaner;
47 
48  if (do_seed_clean)
49  itconf.m_seed_cleaner(seeds, itconf, eoh.refBeamSpot());
50 
51  // Check nans in seeds -- this should not be needed when Slava fixes
52  // the track parameter coordinate transformation.
53  builder.seed_post_cleaning(seeds);
54 
55  if (itconf.m_requires_seed_hit_sorting) {
56  for (auto &s : seeds)
57  s.sortHitsByLayer(); // sort seed hits for the matched hits (I hope it works here)
58  }
59 
60  builder.find_tracks_load_seeds(seeds, do_seed_clean);
61 
62  builder.findTracksCloneEngine();
63 
64  // Pre backward-fit filtering.
65  filter_candidates_func pre_filter;
66  if (do_backward_fit && itconf.m_pre_bkfit_filter)
67  pre_filter = [&](const TrackCand &tc, const MkJob &jb) -> bool {
68  return itconf.m_pre_bkfit_filter(tc, jb) && StdSeq::qfilter_nan_n_silly<TrackCand>(tc, jb);
69  };
70  else if (itconf.m_pre_bkfit_filter)
71  pre_filter = itconf.m_pre_bkfit_filter;
72  else if (do_backward_fit)
73  pre_filter = StdSeq::qfilter_nan_n_silly<TrackCand>;
74  // pre_filter can be null if we are not doing backward fit as nan_n_silly will be run below.
75  if (pre_filter)
76  builder.filter_comb_cands(pre_filter, true);
77 
78  job.switch_to_backward();
79 
80  if (do_backward_fit) {
81  if (itconf.m_backward_search) {
82  builder.compactifyHitStorageForBestCand(itconf.m_backward_drop_seed_hits, itconf.m_backward_fit_min_hits);
83  }
84 
85  builder.backwardFit();
86 
87  if (itconf.m_backward_search) {
88  builder.beginBkwSearch();
89  builder.findTracksCloneEngine(SteeringParams::IT_BkwSearch);
90  }
91  }
92 
93  // Post backward-fit filtering.
94  filter_candidates_func post_filter;
95  if (do_backward_fit && itconf.m_post_bkfit_filter)
96  post_filter = [&](const TrackCand &tc, const MkJob &jb) -> bool {
97  return itconf.m_post_bkfit_filter(tc, jb) && StdSeq::qfilter_nan_n_silly<TrackCand>(tc, jb);
98  };
99  else
100  post_filter = StdSeq::qfilter_nan_n_silly<TrackCand>;
101  // post_filter is always at least doing nan_n_silly filter.
102  builder.filter_comb_cands(post_filter, true);
103 
104  if (do_backward_fit && itconf.m_backward_search)
105  builder.endBkwSearch();
106 
107  builder.export_best_comb_cands(out_tracks, true);
108 
109  if (do_remove_duplicates && itconf.m_duplicate_cleaner) {
110  itconf.m_duplicate_cleaner(out_tracks, itconf);
111  }
112 
113  builder.end_event();
114  builder.release_memory();
115  }
std::function< filter_candidates_cf > filter_candidates_func
Definition: FunctionTypes.h:28

◆ runBtpCe_MultiIter()

std::vector< double > mkfit::runBtpCe_MultiIter ( Event ev,
const EventOfHits eoh,
MkBuilder builder,
int  n 
)

Definition at line 417 of file buildtestMPlex.cc.

References analysisFilters_cff::algorithms, mkfit::Config::backwardFit, mkfit::MkBuilder::backwardFit(), mkfit::Config::backwardSearch, mkfit::MkBuilder::begin_event(), mkfit::MkBuilder::beginBkwSearch(), nano_mu_local_reco_cff::bool, mkfit::Config::cmssw_val, mkfit::MkBuilder::compactifyHitStorageForBestCand(), dtime(), mkfit::MkBuilder::end_event(), mkfit::MkBuilder::endBkwSearch(), makeMEIFBenchmarkPlots::ev, mkfit::MkBuilder::export_tracks(), mkfit::MkBuilder::filter_comb_cands(), spr::find(), mkfit::MkBuilder::find_tracks_load_seeds(), mkfit::MkBuilder::findTracksCloneEngine(), mkfit::SteeringParams::IT_BkwSearch, mkfit::Config::ItrInfo, mkfit::IterationConfig::m_backward_drop_seed_hits, mkfit::IterationConfig::m_backward_fit_min_hits, mkfit::IterationConfig::m_backward_search, mkfit::IterationConfig::m_duplicate_cleaner, mkfit::IterationMaskIfc::m_mask_vector, mkfit::IterationConfig::m_post_bkfit_filter, mkfit::IterationConfig::m_pre_bkfit_filter, mkfit::IterationConfig::m_requires_seed_hit_sorting, mkfit::IterationConfig::m_seed_cleaner, mkfit::IterationConfig::m_track_algorithm, eostools::move(), dqmiodumpmetadata::n, mkfit::StdSeq::prep_simtracks(), mkfit::StdSeq::Quality::quality_val(), mkfit::Config::quality_val, mkfit::MkBuilder::ref_tracks_nc(), mkfit::EventOfHits::refBeamSpot(), mkfit::MkBuilder::release_memory(), mkfit::StdSeq::root_val(), alignCSCRings::s, mkfit::MkBuilder::seed_post_cleaning(), DetachedQuadStep_cff::seeds, mkfit::MkBuilder::select_best_comb_cands(), mkfit::Config::sim_val, submitPVValidationJobs::t, hcalRecHitTable_cff::time, and mkfit::Config::TrkInfo.

Referenced by test_standard().

417  {
418  std::vector<double> timevec;
419  if (n <= 0)
420  return timevec;
421  timevec.resize(n + 1, 0.0);
422 
423  const bool validation_on = (Config::sim_val || Config::quality_val);
424 
425  TrackVec seeds_used;
426  TrackVec seeds1;
427 
428  if (validation_on) {
429  for (auto const &s : ev.seedTracks_) {
430  //keep seeds form the first n iterations for processing
431  if (std::find(algorithms, algorithms + n, s.algoint()) != algorithms + n)
432  seeds1.push_back(s);
433  }
434  ev.seedTracks_.swap(seeds1); //necessary for the validation - PrepareSeeds
435  ev.relabel_bad_seedtracks(); //necessary for the validation - PrepareSeeds
436  }
437 
438  IterationMaskIfc mask_ifc;
439  TrackVec seeds;
440  TrackVec tmp_tvec;
441 
442  for (int it = 0; it <= n - 1; ++it) {
443  const IterationConfig &itconf = Config::ItrInfo[it];
444 
445  // To disable hit-masks, pass nullptr in place of &mask_ifc to MkJob ctor
446  // and optionally comment out ev.fill_hitmask_bool_vectors() call.
447 
448  ev.fill_hitmask_bool_vectors(itconf.m_track_algorithm, mask_ifc.m_mask_vector);
449 
450  MkJob job({Config::TrkInfo, itconf, eoh, eoh.refBeamSpot(), &mask_ifc});
451 
452  builder.begin_event(&job, &ev, __func__);
453 
454  { // We could partition seeds once, store beg, end for each iteration in a map or vector.
455  seeds.clear();
456  int nc = 0;
457  for (auto &s : ev.seedTracks_) {
458  if (s.algoint() == itconf.m_track_algorithm) {
459  if (itconf.m_requires_seed_hit_sorting) {
460  s.sortHitsByLayer();
461  }
462  seeds.push_back(s);
463  ++nc;
464  } else if (nc > 0)
465  break;
466  }
467  }
468 
469  bool do_seed_clean = bool(itconf.m_seed_cleaner);
470 
471  if (do_seed_clean)
472  itconf.m_seed_cleaner(seeds, itconf, eoh.refBeamSpot());
473 
474  builder.seed_post_cleaning(seeds);
475 
476  // Add protection in case no seeds are found for iteration
477  if (seeds.size() <= 0)
478  continue;
479  ev.setCurrentSeedTracks(seeds);
480  ev.simLabelForCurrentSeed(0);
481 
482  builder.find_tracks_load_seeds(seeds, do_seed_clean);
483 
484  double time = dtime();
485 
486  builder.findTracksCloneEngine();
487 
488  timevec[it] = dtime() - time;
489  timevec[n] += timevec[it];
490 
491  // Print min and max size of hots vectors of CombCands.
492  // builder.find_min_max_hots_size();
493 
494  if (validation_on)
495  seeds_used.insert(seeds_used.end(), seeds.begin(), seeds.end()); //cleaned seeds need to be stored somehow
496 
497  // Pre backward-fit filtering.
498  // Note -- slightly different logic than run_OneIteration as we always do nan filters for
499  // export for validation.
500  filter_candidates_func pre_filter;
501  if (itconf.m_pre_bkfit_filter)
502  pre_filter = [&](const TrackCand &tc, const MkJob &jb) -> bool {
503  return itconf.m_pre_bkfit_filter(tc, jb) && StdSeq::qfilter_nan_n_silly<TrackCand>(tc, jb);
504  };
505  else
506  pre_filter = StdSeq::qfilter_nan_n_silly<TrackCand>;
507  // pre_filter is always at least doing nan_n_silly filter.
508  builder.filter_comb_cands(pre_filter, true);
509 
510  builder.select_best_comb_cands();
511 
512  {
513  builder.export_tracks(tmp_tvec);
514  if (itconf.m_duplicate_cleaner)
515  itconf.m_duplicate_cleaner(builder.ref_tracks_nc(), itconf);
516  ev.candidateTracks_.reserve(ev.candidateTracks_.size() + tmp_tvec.size());
517  for (auto &&t : tmp_tvec)
518  ev.candidateTracks_.emplace_back(std::move(t));
519  tmp_tvec.clear();
520  }
521 
522  job.switch_to_backward();
523 
524  // now do backwards fit... do we want to time this section?
525  if (Config::backwardFit) {
526  // a) TrackVec version:
527  // builder.backwardFitBH();
528 
529  // b) Version that runs on CombCand / TrackCand
530  const bool do_backward_search = Config::backwardSearch && itconf.m_backward_search;
531 
532  // We copy seed-hits into Candidates ... now we have to remove them so backward fit stops
533  // before reaching seeding region. Ideally, we wouldn't add them in the first place but
534  // if we want to export full tracks above we need to hold on to them (alternatively, we could
535  // have a pointer to seed track in CombCandidate and copy them from there).
536  if (do_backward_search)
537  builder.compactifyHitStorageForBestCand(itconf.m_backward_drop_seed_hits, itconf.m_backward_fit_min_hits);
538 
539  builder.backwardFit();
540 
541  if (do_backward_search) {
542  builder.beginBkwSearch();
543  builder.findTracksCloneEngine(SteeringParams::IT_BkwSearch);
544  }
545 
546  // Post backward-fit filtering.
547  // Note -- slightly different logic than run_OneIteration as we export both pre and post
548  // backward-fit tracks.
549  filter_candidates_func post_filter;
550  if (itconf.m_post_bkfit_filter)
551  post_filter = [&](const TrackCand &tc, const MkJob &jb) -> bool {
552  return itconf.m_post_bkfit_filter(tc, jb) && StdSeq::qfilter_nan_n_silly<TrackCand>(tc, jb);
553  };
554  else
555  post_filter = StdSeq::qfilter_nan_n_silly<TrackCand>;
556  // post_filter is always at least doing nan_n_silly filter.
557  builder.filter_comb_cands(post_filter, true);
558 
559  if (do_backward_search)
560  builder.endBkwSearch();
561 
562  builder.select_best_comb_cands(true); // true -> clear m_tracks as they were already filled once above
563 
564  if (itconf.m_duplicate_cleaner)
565  itconf.m_duplicate_cleaner(builder.ref_tracks_nc(), itconf);
566 
567  builder.export_tracks(ev.fitTracks_);
568  }
569  ev.resetCurrentSeedTracks();
570 
571  builder.end_event();
572  }
573 
574  // MIMI - Fake back event pointer for final processing (that should be done elsewhere)
575  MkJob job({Config::TrkInfo, Config::ItrInfo[0], eoh, eoh.refBeamSpot()});
576  builder.begin_event(&job, &ev, __func__);
577 
578  if (validation_on) {
580  //swap for the cleaned seeds
581  ev.seedTracks_.swap(seeds_used);
582  }
583 
584  check_nan_n_silly_candidates(ev);
585 
587  check_nan_n_silly_bkfit(ev);
588 
589  // validation section
590  if (Config::quality_val) {
591  StdSeq::Quality qval;
592  qval.quality_val(&ev);
593  } else if (Config::sim_val || Config::cmssw_val) {
595  }
596 
597  // ev.print_tracks(ev.candidateTracks_, true);
598 
599  // MIMI Unfake.
600  builder.end_event();
601 
602  // In CMSSW runOneIter we now release memory for comb-cands:
603  builder.release_memory();
604 
605  return timevec;
606  }
void root_val(Event *event)
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:19
void prep_simtracks(Event *event)
std::function< filter_candidates_cf > filter_candidates_func
Definition: FunctionTypes.h:28
TrackerInfo TrkInfo
IterationsInfo ItrInfo
std::vector< Track > TrackVec
double dtime()
def move(src, dest)
Definition: eostools.py:511

◆ runBuildingTestPlexBestHit()

double mkfit::runBuildingTestPlexBestHit ( Event ev,
const EventOfHits eoh,
MkBuilder builder 
)

Definition at line 108 of file buildtestMPlex.cc.

References analysisFilters_cff::algorithms, mkfit::Config::backwardFit, mkfit::MkBuilder::backwardFitBH(), mkfit::MkBuilder::begin_event(), mkfit::Config::cmssw_val, dtime(), mkfit::MkBuilder::end_event(), makeMEIFBenchmarkPlots::ev, spr::find(), mkfit::MkBuilder::find_tracks_load_seeds_BH(), mkfit::MkBuilder::findTracksBestHit(), mkfit::Config::ItrInfo, mkfit::IterationMaskIfc::m_mask_vector, mkfit::IterationConfig::m_track_algorithm, mkfit::StdSeq::Quality::quality_val(), mkfit::Config::quality_val, mkfit::MkBuilder::ref_tracks(), mkfit::EventOfHits::refBeamSpot(), mkfit::StdSeq::root_val(), alignCSCRings::s, mkfit::Config::sim_val, hcalRecHitTable_cff::time, and mkfit::Config::TrkInfo.

Referenced by test_standard().

108  {
109  const IterationConfig &itconf = Config::ItrInfo[0];
110 
111  const bool validation_on = (Config::sim_val || Config::quality_val);
112 
113  TrackVec seeds1;
114  if (validation_on) {
115  unsigned int algorithms[] = {4}; //only initialStep
116 
117  for (auto const &s : ev.seedTracks_) {
118  //keep seeds form the first iteration for processing
119  if (std::find(algorithms, algorithms + 1, s.algoint()) != algorithms + 1)
120  seeds1.push_back(s);
121  }
122  ev.seedTracks_.swap(seeds1); //necessary for the validation - PrepareSeeds
123  ev.relabel_bad_seedtracks(); //necessary for the validation - PrepareSeeds
124  }
125 
126  IterationMaskIfc mask_ifc;
127 
128  // To disable hit-masks, pass nullptr in place of &mask_ifc to MkJob ctor
129  // and optionally comment out ev.fill_hitmask_bool_vectors() call.
130 
131  ev.fill_hitmask_bool_vectors(itconf.m_track_algorithm, mask_ifc.m_mask_vector);
132 
133  MkJob job({Config::TrkInfo, itconf, eoh, eoh.refBeamSpot(), &mask_ifc});
134 
135  builder.begin_event(&job, &ev, __func__);
136 
137  bool seeds_sorted = false;
138  // CCCC builder.PrepareSeeds();
139 
140  // EventOfCandidates event_of_cands;
141  builder.find_tracks_load_seeds_BH(ev.seedTracks_, seeds_sorted);
142 
143 #ifdef USE_VTUNE_PAUSE
144  __SSC_MARK(0x111); // use this to resume Intel SDE at the same point
145  __itt_resume();
146 #endif
147 
148  double time = dtime();
149 
150  builder.findTracksBestHit();
151 
152  time = dtime() - time;
153 
154 #ifdef USE_VTUNE_PAUSE
155  __itt_pause();
156  __SSC_MARK(0x222); // use this to pause Intel SDE at the same point
157 #endif
158 
159  // Hack, get the tracks out.
160  ev.candidateTracks_ = builder.ref_tracks();
161 
162  // For best hit, the candidateTracks_ vector is the direct input to the backward fit so only need to do clean_duplicates once
164  //Mark tracks as duplicates; if within CMSSW, remove duplicate tracks before backward fit
165  // CCCC if (Config::removeDuplicates) {
166  // CCCC StdSeq::clean_duplicates(ev.candidateTracks_);
167  // CCCC }
168  }
169 
170  job.switch_to_backward();
171 
172  // now do backwards fit... do we want to time this section?
173  if (Config::backwardFit) {
174  builder.backwardFitBH();
175  ev.fitTracks_ = builder.ref_tracks();
176  }
177 
178  if (Config::quality_val) {
179  StdSeq::Quality qval;
180  qval.quality_val(&ev);
181  } else if (Config::sim_val || Config::cmssw_val) {
183  }
184 
185  builder.end_event();
186 
187  // ev.print_tracks(ev.candidateTracks_, true);
188 
189  if (validation_on) {
190  ev.seedTracks_.swap(seeds1);
191  }
192 
193  return time;
194  }
void root_val(Event *event)
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:19
TrackerInfo TrkInfo
IterationsInfo ItrInfo
std::vector< Track > TrackVec
double dtime()

◆ runBuildingTestPlexCloneEngine()

double mkfit::runBuildingTestPlexCloneEngine ( Event ev,
const EventOfHits eoh,
MkBuilder builder 
)

Definition at line 297 of file buildtestMPlex.cc.

References analysisFilters_cff::algorithms, mkfit::Config::backwardFit, mkfit::MkBuilder::backwardFitBH(), mkfit::MkBuilder::begin_event(), mkfit::Config::cmssw_val, dtime(), mkfit::MkBuilder::end_event(), makeMEIFBenchmarkPlots::ev, mkfit::MkBuilder::export_best_comb_cands(), spr::find(), mkfit::MkBuilder::find_tracks_load_seeds(), mkfit::MkBuilder::findTracksCloneEngine(), mkfit::Config::ItrInfo, mkfit::IterationMaskIfc::m_mask_vector, mkfit::IterationConfig::m_track_algorithm, mkfit::StdSeq::Quality::quality_val(), mkfit::Config::quality_val, mkfit::MkBuilder::ref_tracks(), mkfit::EventOfHits::refBeamSpot(), mkfit::StdSeq::root_val(), alignCSCRings::s, mkfit::MkBuilder::select_best_comb_cands(), mkfit::Config::sim_val, hcalRecHitTable_cff::time, and mkfit::Config::TrkInfo.

Referenced by test_standard().

297  {
298  const IterationConfig &itconf = Config::ItrInfo[0];
299 
300  const bool validation_on = (Config::sim_val || Config::quality_val);
301 
302  TrackVec seeds1;
303  if (validation_on) {
304  unsigned int algorithms[] = {4}; //only initialStep
305 
306  for (auto const &s : ev.seedTracks_) {
307  //keep seeds form the first iteration for processing
308  if (std::find(algorithms, algorithms + 1, s.algoint()) != algorithms + 1)
309  seeds1.push_back(s);
310  }
311  ev.seedTracks_.swap(seeds1); //necessary for the validation - PrepareSeeds
312  ev.relabel_bad_seedtracks(); //necessary for the validation - PrepareSeeds
313  }
314 
315  IterationMaskIfc mask_ifc;
316 
317  // To disable hit-masks, pass nullptr in place of &mask_ifc to MkJob ctor
318  // and optionally comment out ev.fill_hitmask_bool_vectors() call.
319 
320  ev.fill_hitmask_bool_vectors(itconf.m_track_algorithm, mask_ifc.m_mask_vector);
321 
322  MkJob job({Config::TrkInfo, itconf, eoh, eoh.refBeamSpot(), &mask_ifc});
323 
324  builder.begin_event(&job, &ev, __func__);
325 
326  bool seeds_sorted = false;
327  // CCCC builder.PrepareSeeds();
328  ev.setCurrentSeedTracks(ev.seedTracks_);
329 
330  builder.find_tracks_load_seeds(ev.seedTracks_, seeds_sorted);
331 
332 #ifdef USE_VTUNE_PAUSE
333  __SSC_MARK(0x111); // use this to resume Intel SDE at the same point
334  __itt_resume();
335 #endif
336 
337  double time = dtime();
338 
339  builder.findTracksCloneEngine();
340 
341  time = dtime() - time;
342 
343 #ifdef USE_VTUNE_PAUSE
344  __itt_pause();
345  __SSC_MARK(0x222); // use this to pause Intel SDE at the same point
346 #endif
347 
348  check_nan_n_silly_candidates(ev);
349 
350  // first store candidate tracks - needed for BH backward fit and root_validation
351  ev.candidateTracks_.clear();
352  builder.export_best_comb_cands(ev.candidateTracks_);
353 
354  job.switch_to_backward();
355 
356  // now do backwards fit... do we want to time this section?
357  if (Config::backwardFit) {
358  // a) TrackVec version:
359  builder.select_best_comb_cands();
360  builder.backwardFitBH();
361  ev.fitTracks_ = builder.ref_tracks();
362 
363  // b) Version that runs on CombCand / TrackCand
364  // builder.backwardFit();
365  // builder.quality_store_tracks(ev.fitTracks_);
366 
367  check_nan_n_silly_bkfit(ev);
368  }
369 
370  // CCCC StdSeq::handle_duplicates(&ev);
371 
372  // validation section
373  if (Config::quality_val) {
374  StdSeq::Quality qval;
375  qval.quality_val(&ev);
376  } else if (Config::sim_val || Config::cmssw_val) {
378  }
379 
380  ev.resetCurrentSeedTracks();
381 
382  builder.end_event();
383 
384  // ev.print_tracks(ev.candidateTracks_, true);
385 
386  if (validation_on) {
387  ev.seedTracks_.swap(seeds1);
388  }
389 
390  return time;
391  }
void root_val(Event *event)
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:19
TrackerInfo TrkInfo
IterationsInfo ItrInfo
std::vector< Track > TrackVec
double dtime()

◆ runBuildingTestPlexDumbCMSSW()

void mkfit::runBuildingTestPlexDumbCMSSW ( Event ev,
const EventOfHits eoh,
MkBuilder builder 
)

Definition at line 90 of file buildtestMPlex.cc.

References mkfit::MkBuilder::begin_event(), mkfit::MkBuilder::end_event(), makeMEIFBenchmarkPlots::ev, mkfit::Config::ItrInfo, mkfit::EventOfHits::refBeamSpot(), mkfit::StdSeq::root_val_dumb_cmssw(), mkfit::Config::sim_val_for_cmssw, and mkfit::Config::TrkInfo.

Referenced by test_standard().

90  {
91  const IterationConfig &itconf = Config::ItrInfo[0];
92 
93  MkJob job({Config::TrkInfo, itconf, eoh, eoh.refBeamSpot()});
94 
95  builder.begin_event(&job, &ev, __func__);
96 
99  }
100 
101  builder.end_event();
102  }
TrackerInfo TrkInfo
IterationsInfo ItrInfo
void root_val_dumb_cmssw(Event *event)

◆ runBuildingTestPlexStandard()

double mkfit::runBuildingTestPlexStandard ( Event ev,
const EventOfHits eoh,
MkBuilder builder 
)

Definition at line 200 of file buildtestMPlex.cc.

References analysisFilters_cff::algorithms, mkfit::Config::backwardFit, mkfit::MkBuilder::backwardFitBH(), mkfit::MkBuilder::begin_event(), mkfit::Config::cmssw_val, dtime(), mkfit::MkBuilder::end_event(), makeMEIFBenchmarkPlots::ev, mkfit::MkBuilder::export_best_comb_cands(), spr::find(), mkfit::MkBuilder::find_tracks_load_seeds(), mkfit::MkBuilder::findTracksStandard(), mkfit::Config::ItrInfo, mkfit::IterationMaskIfc::m_mask_vector, mkfit::IterationConfig::m_track_algorithm, mkfit::StdSeq::Quality::quality_val(), mkfit::Config::quality_val, mkfit::MkBuilder::ref_tracks(), mkfit::EventOfHits::refBeamSpot(), mkfit::StdSeq::root_val(), alignCSCRings::s, mkfit::MkBuilder::select_best_comb_cands(), mkfit::Config::sim_val, hcalRecHitTable_cff::time, and mkfit::Config::TrkInfo.

Referenced by test_standard().

200  {
201  const IterationConfig &itconf = Config::ItrInfo[0];
202 
203  const bool validation_on = (Config::sim_val || Config::quality_val);
204 
205  TrackVec seeds1;
206  if (validation_on) {
207  unsigned int algorithms[] = {4}; //only initialStep
208 
209  for (auto const &s : ev.seedTracks_) {
210  //keep seeds form the first iteration for processing
211  if (std::find(algorithms, algorithms + 1, s.algoint()) != algorithms + 1)
212  seeds1.push_back(s);
213  }
214  ev.seedTracks_.swap(seeds1); //necessary for the validation - PrepareSeeds
215  ev.relabel_bad_seedtracks(); //necessary for the validation - PrepareSeeds
216  }
217 
218  IterationMaskIfc mask_ifc;
219 
220  // To disable hit-masks, pass nullptr in place of &mask_ifc to MkJob ctor
221  // and optionally comment out ev.fill_hitmask_bool_vectors() call.
222 
223  ev.fill_hitmask_bool_vectors(itconf.m_track_algorithm, mask_ifc.m_mask_vector);
224 
225  MkJob job({Config::TrkInfo, itconf, eoh, eoh.refBeamSpot(), &mask_ifc});
226 
227  builder.begin_event(&job, &ev, __func__);
228 
229  bool seeds_sorted = false;
230  // CCCC builder.PrepareSeeds();
231  ev.setCurrentSeedTracks(ev.seedTracks_);
232  ev.simLabelForCurrentSeed(0);
233 
234  builder.find_tracks_load_seeds(ev.seedTracks_, seeds_sorted);
235 
236 #ifdef USE_VTUNE_PAUSE
237  __SSC_MARK(0x111); // use this to resume Intel SDE at the same point
238  __itt_resume();
239 #endif
240 
241  double time = dtime();
242 
243  builder.findTracksStandard();
244 
245  time = dtime() - time;
246 
247 #ifdef USE_VTUNE_PAUSE
248  __itt_pause();
249  __SSC_MARK(0x222); // use this to pause Intel SDE at the same point
250 #endif
251 
252  check_nan_n_silly_candidates(ev);
253 
254  // first store candidate tracks
255  ev.candidateTracks_.clear();
256  builder.export_best_comb_cands(ev.candidateTracks_);
257 
258  job.switch_to_backward();
259 
260  // now do backwards fit... do we want to time this section?
261  if (Config::backwardFit) {
262  // Using the TrackVec version until we home in on THE backward fit etc.
263  // builder.backwardFit();
264  builder.select_best_comb_cands();
265  builder.backwardFitBH();
266  ev.fitTracks_ = builder.ref_tracks();
267 
268  check_nan_n_silly_bkfit(ev);
269  }
270 
271  // CCCC StdSeq::handle_duplicates(&ev);
272 
273  if (Config::quality_val) {
274  StdSeq::Quality qval;
275  qval.quality_val(&ev);
276  } else if (Config::sim_val || Config::cmssw_val) {
278  }
279 
280  ev.resetCurrentSeedTracks();
281 
282  builder.end_event();
283 
284  // ev.print_tracks(ev.candidateTracks_, true);
285 
286  if (validation_on) {
287  ev.seedTracks_.swap(seeds1);
288  }
289 
290  return time;
291  }
void root_val(Event *event)
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:19
TrackerInfo TrkInfo
IterationsInfo ItrInfo
std::vector< Track > TrackVec
double dtime()

◆ runFittingTestPlex()

double mkfit::runFittingTestPlex ( Event ev,
std::vector< Track > &  rectracks 
)

Definition at line 40 of file fittestMPlex.cc.

References submitPVResolutionJobs::count, dtime(), makeMEIFBenchmarkPlots::ev, mkfit::Config::fit_val, g_exe_ctx, mps_fire::i, mkfit::ExecutionContext::m_fitters, SiStripPI::max, mkfit::Config::nLayers, NN, mkfit::Config::numSeedsPerTask, mkfit::Config::numThreadsFinder, mkfit::ExecutionContext::populate(), and hcalRecHitTable_cff::time.

40  {
42  std::vector<Track>& simtracks = ev.simTracks_;
43 
44  const int Nhits = Config::nLayers;
45  // XXX What if there's a missing / double layer?
46  // Eventually, should sort track vector by number of hits!
47  // And pass the number in on each "setup" call.
48  // Reserves should be made for maximum possible number (but this is just
49  // measurments errors, params).
50 
51  int theEnd = simtracks.size();
52  int count = (theEnd + NN - 1) / NN;
53 
54 #ifdef USE_VTUNE_PAUSE
55  __SSC_MARK(0x111); // use this to resume Intel SDE at the same point
56  __itt_resume();
57 #endif
58 
59  double time = dtime();
60 
61  tbb::parallel_for(tbb::blocked_range<int>(0, count, std::max(1, Config::numSeedsPerTask / NN)),
62  [&](const tbb::blocked_range<int>& i) {
63  std::unique_ptr<MkFitter, decltype(retfitr)> mkfp(g_exe_ctx.m_fitters.GetFromPool(), retfitr);
64  mkfp->setNhits(Nhits);
65  for (int it = i.begin(); it < i.end(); ++it) {
66  int itrack = it * NN;
67  int end = itrack + NN;
68  /*
69  * MT, trying to slurp and fit at the same time ...
70  if (theEnd < end) {
71  end = theEnd;
72  mkfp->inputTracksAndHits(simtracks, ev.layerHits_, itrack, end);
73  } else {
74  mkfp->slurpInTracksAndHits(simtracks, ev.layerHits_, itrack, end); // only safe for a full matriplex
75  }
76 
77  if (Config::cf_fitting) mkfp->ConformalFitTracks(true, itrack, end);
78  mkfp->FitTracks(end - itrack, &ev, true);
79  */
80 
81  mkfp->inputTracksForFit(simtracks, itrack, end);
82 
83  // XXXX MT - for this need 3 points in ... right
84  // XXXX if (Config::cf_fitting) mkfp->ConformalFitTracks(true, itrack, end);
85 
86  mkfp->fitTracksWithInterSlurp(ev.layerHits_, end - itrack);
87 
88  mkfp->outputFittedTracks(rectracks, itrack, end);
89  }
90  });
91 
92  time = dtime() - time;
93 
94 #ifdef USE_VTUNE_PAUSE
95  __itt_pause();
96  __SSC_MARK(0x222); // use this to pause Intel SDE at the same point
97 #endif
98 
99  if (Config::fit_val)
100  ev.validate();
101 
102  return time;
103  }
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
ExecutionContext g_exe_ctx
Definition: MkBuilder.cc:49
constexpr int numThreadsFinder
Definition: Config.h:87
Pool< MkFitter > m_fitters
Definition: MkBuilder.cc:39
void populate(int n_thr)
Definition: MkBuilder.cc:42
double dtime()
constexpr int numSeedsPerTask
Definition: Config.h:89

◆ sincos4()

void mkfit::sincos4 ( const float  x,
float &  sin,
float &  cos 
)
inline

Definition at line 13 of file Matrix.h.

References funct::cos(), funct::sin(), and x.

Referenced by helixAtRFromIterativeCCSFullJac(), and helixAtZ().

13  {
14  // Had this writen with explicit division by factorial.
15  // The *whole* fitting test ran like 2.5% slower on MIC, sigh.
16 
17  const float x2 = x * x;
18  cos = 1.f - 0.5f * x2 + 0.04166667f * x2 * x2;
19  sin = x - 0.16666667f * x * x2;
20  }
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
float x

◆ sortByEta()

bool mkfit::sortByEta ( const Hit hit1,
const Hit hit2 
)
inline

Definition at line 25 of file buildtestMPlex.cc.

References mkfit::Hit::eta().

25 { return hit1.eta() < hit2.eta(); }

◆ sortByHitsChi2() [1/2]

bool mkfit::sortByHitsChi2 ( const std::pair< Track, TrackState > &  cand1,
const std::pair< Track, TrackState > &  cand2 
)
inline

Definition at line 14 of file buildtestMPlex.cc.

14  {
15  if (cand1.first.nFoundHits() == cand2.first.nFoundHits())
16  return cand1.first.chi2() < cand2.first.chi2();
17 
18  return cand1.first.nFoundHits() > cand2.first.nFoundHits();
19  }

◆ sortByHitsChi2() [2/2]

bool mkfit::sortByHitsChi2 ( const Track cand1,
const Track cand2 
)
inline

Definition at line 599 of file Track.h.

References mkfit::TrackBase::chi2(), and mkfit::Track::nFoundHits().

599  {
600  if (cand1.nFoundHits() == cand2.nFoundHits())
601  return cand1.chi2() < cand2.chi2();
602  return cand1.nFoundHits() > cand2.nFoundHits();
603  }

◆ sortByPhi()

bool mkfit::sortByPhi ( const Hit hit1,
const Hit hit2 
)
inline

Definition at line 21 of file buildtestMPlex.cc.

References mkfit::Hit::x(), and mkfit::Hit::y().

21  {
22  return std::atan2(hit1.y(), hit1.x()) < std::atan2(hit2.y(), hit2.x());
23  }

◆ sortByScoreCand()

bool mkfit::sortByScoreCand ( const Track cand1,
const Track cand2 
)
inline

Definition at line 605 of file Track.h.

References mkfit::TrackBase::score().

Referenced by mkfit::TTreeValidation::mapRefTkToRecoTks().

605 { return cand1.score() > cand2.score(); }

◆ sortByScoreStruct()

bool mkfit::sortByScoreStruct ( const IdxChi2List cand1,
const IdxChi2List cand2 
)
inline

Definition at line 607 of file Track.h.

References mkfit::IdxChi2List::score.

607  {
608  return cand1.score > cand2.score;
609  }

◆ sortByScoreTrackCand()

bool mkfit::sortByScoreTrackCand ( const TrackCand cand1,
const TrackCand cand2 
)
inline

Definition at line 255 of file TrackStructures.h.

References mkfit::TrackBase::score().

Referenced by mkfit::CombCandidate::mergeCandsAndBestShortOne(), and mkfit::CandCloner::processSeedRange().

255  {
256  return cand1.score() > cand2.score();
257  }

◆ sortByZ()

bool mkfit::sortByZ ( const Hit hit1,
const Hit hit2 
)
inline

Definition at line 43 of file buildtestMPlex.cc.

References mkfit::Hit::z().

43 { return hit1.z() < hit2.z(); }
double z
global z - AlignmentGeometry::z0, mm
Definition: HitCollection.h:27

◆ sortIDsByChi2()

bool mkfit::sortIDsByChi2 ( const idchi2Pair cand1,
const idchi2Pair cand2 
)
inline

Definition at line 235 of file TrackExtra.cc.

Referenced by mkfit::TrackExtra::setCMSSWTrackIDInfoByTrkParams().

235 { return cand1.second < cand2.second; }

◆ sortTracksByEta()

bool mkfit::sortTracksByEta ( const Track track1,
const Track track2 
)
inline

Definition at line 27 of file buildtestMPlex.cc.

References mkfit::TrackBase::momEta().

27 { return track1.momEta() < track2.momEta(); }

◆ sortTracksByPhi()

bool mkfit::sortTracksByPhi ( const Track track1,
const Track track2 
)
inline

Definition at line 29 of file buildtestMPlex.cc.

References mkfit::TrackBase::momPhi().

29 { return track1.momPhi() < track2.momPhi(); }

◆ sqr()

template<typename T >
T mkfit::sqr ( T  x)
inline

Definition at line 14 of file Hit.h.

References x.

14  {
15  return x * x;
16  }
float x

◆ squashPhiGeneral() [1/2]

float mkfit::squashPhiGeneral ( float  phi)
inline

◆ squashPhiGeneral() [2/2]

template<typename Vector >
void mkfit::squashPhiGeneral ( Vector v)
inline

Definition at line 645 of file Track.h.

References mps_fire::i, squashPhiGeneral(), and findQualityFiles::v.

645  {
646  const int i = v.kSize - 2; // phi index
647  v[i] = squashPhiGeneral(v[i]);
648  }
void squashPhiGeneral(Vector &v)
Definition: Track.h:645

◆ squashPhiMinimal()

float mkfit::squashPhiMinimal ( float  phi)
inline

◆ squashPhiMPlex()

void mkfit::squashPhiMPlex ( MPlexLV par,
const int  N_proc 
)
inline

Definition at line 10 of file PropagationMPlex.h.

References dqmiodumpmetadata::n, NN, mkfit::Const::PI, and mkfit::Const::TwoPI.

Referenced by kalmanOperation(), kalmanOperationEndcap(), kalmanOperationPlane(), propagateHelixToPlaneMPlex(), propagateHelixToRMPlex(), and propagateHelixToZMPlex().

10  {
11 #pragma omp simd
12  for (int n = 0; n < NN; ++n) {
13  if (par(n, 4, 0) >= Const::PI)
14  par(n, 4, 0) -= Const::TwoPI;
15  if (par(n, 4, 0) < -Const::PI)
16  par(n, 4, 0) += Const::TwoPI;
17  }
18  }
constexpr float TwoPI
Definition: Config.h:8
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
#define PI
Definition: QcdUeDQM.h:37

◆ squashPhiMPlexGeneral()

void mkfit::squashPhiMPlexGeneral ( MPlexLV par,
const int  N_proc 
)
inline

Definition at line 20 of file PropagationMPlex.h.

References f, mkfit::Const::InvPI, dqmiodumpmetadata::n, NN, mkfit::Const::PI, and mkfit::Const::TwoPI.

20  {
21 #pragma omp simd
22  for (int n = 0; n < NN; ++n) {
23  par(n, 4, 0) -= std::floor(0.5f * Const::InvPI * (par(n, 4, 0) + Const::PI)) * Const::TwoPI;
24  }
25  }
constexpr float TwoPI
Definition: Config.h:8
constexpr Matriplex::idx_t NN
Definition: Matrix.h:43
double f[11][100]
#define PI
Definition: QcdUeDQM.h:37
constexpr float InvPI
Definition: Config.h:12

◆ to_json() [1/12]

void mkfit::to_json ( nlohmann::ordered_json &  nlohmann_json_j,
const mkfit::LayerControl nlohmann_json_t 
)
inline

Definition at line 37 of file IterationConfig.cc.

◆ to_json() [2/12]

void mkfit::to_json ( nlohmann::json nlohmann_json_j,
const mkfit::LayerControl nlohmann_json_t 
)
inline

Definition at line 37 of file IterationConfig.cc.

Referenced by mkfit::ConfigJson::save_Iterations().

◆ to_json() [3/12]

void mkfit::to_json ( nlohmann::ordered_json &  nlohmann_json_j,
const mkfit::SteeringParams nlohmann_json_t 
)
inline

Definition at line 45 of file IterationConfig.cc.

◆ to_json() [4/12]

void mkfit::to_json ( nlohmann::json nlohmann_json_j,
const mkfit::SteeringParams nlohmann_json_t 
)
inline

Definition at line 45 of file IterationConfig.cc.

◆ to_json() [5/12]

void mkfit::to_json ( nlohmann::ordered_json &  nlohmann_json_j,
const mkfit::IterationLayerConfig nlohmann_json_t 
)
inline

Definition at line 54 of file IterationConfig.cc.

◆ to_json() [6/12]

void mkfit::to_json ( nlohmann::json nlohmann_json_j,
const mkfit::IterationLayerConfig nlohmann_json_t 
)
inline

Definition at line 54 of file IterationConfig.cc.

◆ to_json() [7/12]

void mkfit::to_json ( nlohmann::json nlohmann_json_j,
const mkfit::IterationParams nlohmann_json_t 
)
inline

Definition at line 68 of file IterationConfig.cc.

◆ to_json() [8/12]

void mkfit::to_json ( nlohmann::ordered_json &  nlohmann_json_j,
const mkfit::IterationParams nlohmann_json_t 
)
inline

Definition at line 68 of file IterationConfig.cc.

◆ to_json() [9/12]

void mkfit::to_json ( nlohmann::ordered_json &  nlohmann_json_j,
const mkfit::IterationConfig nlohmann_json_t 
)
inline

Definition at line 101 of file IterationConfig.cc.

◆ to_json() [10/12]

void mkfit::to_json ( nlohmann::json nlohmann_json_j,
const mkfit::IterationConfig nlohmann_json_t 
)
inline

Definition at line 101 of file IterationConfig.cc.

◆ to_json() [11/12]

void mkfit::to_json ( nlohmann::ordered_json &  nlohmann_json_j,
const mkfit::IterationsInfo nlohmann_json_t 
)
inline

Definition at line 104 of file IterationConfig.cc.

References CMS_SA_ALLOW, and mutex.

110 {

◆ to_json() [12/12]

void mkfit::to_json ( nlohmann::json nlohmann_json_j,
const mkfit::IterationsInfo nlohmann_json_t 
)
inline

Definition at line 104 of file IterationConfig.cc.

110 {

Variable Documentation

◆ g_debug

bool mkfit::g_debug = true

◆ g_exe_ctx

ExecutionContext mkfit::g_exe_ctx

◆ HH

constexpr Matriplex::idx_t mkfit::HH = 3

Definition at line 46 of file Matrix.h.

Referenced by L1TConfigDumper::analyze().

◆ LL

constexpr Matriplex::idx_t mkfit::LL = 6

◆ MPlexHitIdxMax

constexpr int mkfit::MPlexHitIdxMax = 16
static

Definition at line 12 of file MkFitter.h.

◆ NN

constexpr Matriplex::idx_t mkfit::NN = 8