diff --git a/Tracking/Acts/FaserActsKalmanFilter/FaserActsKalmanFilter/CircleFitTrackSeedTool.h b/Tracking/Acts/FaserActsKalmanFilter/FaserActsKalmanFilter/CircleFitTrackSeedTool.h index 51fcd477af77c1f576db1954921d36dc62758cba..706441ee485bc01f18191c25aa06ea6d0bb0e7ea 100644 --- a/Tracking/Acts/FaserActsKalmanFilter/FaserActsKalmanFilter/CircleFitTrackSeedTool.h +++ b/Tracking/Acts/FaserActsKalmanFilter/FaserActsKalmanFilter/CircleFitTrackSeedTool.h @@ -64,6 +64,7 @@ private: std::vector<Acts::Vector3> fakePositions; double c0, c1, cx, cy, r, chi2, momentum, charge; + Acts::Vector3 direction; size_t size; Acts::CurvilinearTrackParameters get_params(double origin, Acts::BoundSymMatrix cov) const; @@ -73,7 +74,7 @@ private: double getY(double z); void fit(); void fakeFit(double B=0.55); - void linearFit(); + void linearFit(const std::vector<Acts::Vector2> &points); double dx, dy; double m_MeV2GeV = 1e-3; double m_sigma_x = 0.8; diff --git a/Tracking/Acts/FaserActsKalmanFilter/FaserActsKalmanFilter/LinearFit.h b/Tracking/Acts/FaserActsKalmanFilter/FaserActsKalmanFilter/LinearFit.h index 4531156216fceaab578d37662e514e9005226b5f..a7259dc0825bd5e471fac3b1fe129e2c6a4e8c43 100644 --- a/Tracking/Acts/FaserActsKalmanFilter/FaserActsKalmanFilter/LinearFit.h +++ b/Tracking/Acts/FaserActsKalmanFilter/FaserActsKalmanFilter/LinearFit.h @@ -7,11 +7,7 @@ namespace LinearFit { using Vector2 = Eigen::Matrix<double, 2, 1>; -std::pair <Vector2, Vector2 > linearFit(const std::vector<Amg::Vector3D> &spacePoints) { - std::vector<Vector2> points {}; - for (auto sp: spacePoints) { - points.push_back({sp.z(), sp.x()}); - } +std::pair <Vector2, Vector2 > linearFit(const std::vector<Vector2> &points) { size_t nPoints = points.size(); Eigen::Matrix< Vector2::Scalar, Eigen::Dynamic, Eigen::Dynamic > centers(nPoints, 2); for (size_t i = 0; i < nPoints; ++i) centers.row(i) = points[i]; diff --git a/Tracking/Acts/FaserActsKalmanFilter/python/CKF2Config.py b/Tracking/Acts/FaserActsKalmanFilter/python/CKF2Config.py index 41550fcfe4789c7b0fbb226a43bb0b3242f3cf45..3f455e60018deb6e157d5d895a44c95203de5703 100644 --- a/Tracking/Acts/FaserActsKalmanFilter/python/CKF2Config.py +++ b/Tracking/Acts/FaserActsKalmanFilter/python/CKF2Config.py @@ -34,6 +34,11 @@ def CKF2Cfg(flags, **kwargs): acts_tracking_geometry_svc = ActsTrackingGeometrySvcCfg(flags) acc.merge(acts_tracking_geometry_svc ) + if flags.GeoModel.FaserVersion is "FASERNU-01": + origin = 0 + else: + origin = -1900 + # track_seed_tool = CompFactory.ClusterTrackSeedTool() # track_seed_tool = CompFactory.ActsTrackSeedTool() # track_seed_tool = CompFactory.MyTrackSeedTool() @@ -54,7 +59,7 @@ def CKF2Cfg(flags, **kwargs): track_seed_tool.covTheta = initial_variance_inflation[3] * sigma_theta * sigma_theta track_seed_tool.covQOverP = initial_variance_inflation[4] * sigma_qop * sigma_qop track_seed_tool.std_cluster = 0.0231 - track_seed_tool.origin = 0 + track_seed_tool.origin = origin track_seed_tool.TrackCollection = "Segments" trajectory_states_writer_tool = CompFactory.RootTrajectoryStatesWriterTool() @@ -89,20 +94,24 @@ def CKF2Cfg(flags, **kwargs): kalman_fitter1.noDiagnostics = kwargs["noDiagnostics"] kalman_fitter1.ActsLogging = "INFO" kalman_fitter1.SummaryWriter = True - kalman_fitter1.StatesWriter = True + kalman_fitter1.StatesWriter = False kalman_fitter1.SeedCovarianceScale = 10 + kalman_fitter1.isMC = flags.Input.isMC kalman_fitter1.RootTrajectoryStatesWriterTool = trajectory_states_writer_tool1 kalman_fitter1.RootTrajectorySummaryWriterTool = trajectory_summary_writer_tool1 + kalman_fitter1.origin = origin ckf.KalmanFitterTool1 = kalman_fitter1 kalman_fitter2 = CompFactory.KalmanFitterTool(name="fitterTool2", **kwargs) kalman_fitter2.noDiagnostics = kwargs["noDiagnostics"] kalman_fitter2.ActsLogging = "INFO" kalman_fitter2.SummaryWriter = True - kalman_fitter2.StatesWriter = True - kalman_fitter2.SeedCovarianceScale = 1 + kalman_fitter2.StatesWriter = False + kalman_fitter2.SeedCovarianceScale = 10 + kalman_fitter2.isMC = flags.Input.isMC kalman_fitter2.RootTrajectoryStatesWriterTool = trajectory_states_writer_tool2 kalman_fitter2.RootTrajectorySummaryWriterTool = trajectory_summary_writer_tool2 + kalman_fitter2.origin = origin ckf.KalmanFitterTool2 = kalman_fitter2 ckf.TrackSeed = track_seed_tool @@ -110,9 +119,14 @@ def CKF2Cfg(flags, **kwargs): ckf.RootTrajectoryStatesWriterTool = trajectory_states_writer_tool ckf.RootTrajectorySummaryWriterTool = trajectory_summary_writer_tool ckf.PerformanceWriterTool = trajectory_performance_writer_tool + ckf.isMC = flags.Input.isMC + ckf.origin = origin + ckf.SummaryWriter = True + ckf.StatesWriter = False + ckf.PerformanceWriter = False ckf.nMax = 10 ckf.chi2Max = 25 acc.addEventAlgo(ckf) - acc.merge(CKF2_OutputCfg(flags)) + # acc.merge(CKF2_OutputCfg(flags)) return acc diff --git a/Tracking/Acts/FaserActsKalmanFilter/src/CircleFitTrackSeedTool.cxx b/Tracking/Acts/FaserActsKalmanFilter/src/CircleFitTrackSeedTool.cxx index 5cc5aa20a80dbab1ac2120c48536b2652ee3ac1b..ab6a54b5e946093752439c16c318872976862d41 100644 --- a/Tracking/Acts/FaserActsKalmanFilter/src/CircleFitTrackSeedTool.cxx +++ b/Tracking/Acts/FaserActsKalmanFilter/src/CircleFitTrackSeedTool.cxx @@ -102,20 +102,23 @@ StatusCode CircleFitTrackSeedTool::run() { std::vector<Segment> combination {}; std::vector<Seed> seeds {}; // create seeds from four stations - // go(segments, combination, seeds, 0, 4); - // create seeds from three stations go(segments, combination, seeds, 0, 4); if (seeds.size() < 2) { + // create seeds from three stations go(segments, combination, seeds, 0, 3); } // create seeds from two stations if (seeds.size() < 2) { go(segments, combination, seeds, 0, 2); } + if (seeds.size() < 2) { + go(segments, combination, seeds, 0, 1); + } std::list<Seed> allSeeds; for (const Seed &seed : seeds) allSeeds.push_back(seed); + /* allSeeds.sort([](const Seed &left, const Seed &right) { if (left.size > right.size) return true; if (left.size < right.size) return false; @@ -131,7 +134,12 @@ StatusCode CircleFitTrackSeedTool::run() { return (p.size <= 12) || ((p.clusterSet & selected.clusterSet).count() > 6); }); } + */ + std::vector<Seed> selectedSeeds {}; + for (const Seed &seed : allSeeds) { + selectedSeeds.push_back(seed); + } Acts::BoundSymMatrix cov = Acts::BoundSymMatrix::Zero(); cov(Acts::eBoundLoc0, Acts::eBoundLoc0) = m_covLoc0; @@ -240,11 +248,32 @@ CircleFitTrackSeedTool::Seed::Seed(const std::vector<Segment> &segments) : clusters.insert(clusters.end(), seg.clusters.begin(), seg.clusters.end()); spacePoints.insert(spacePoints.end(), seg.spacePoints.begin(), seg.spacePoints.end()); positions.push_back(seg.position); + // TODO use reconstruct space points instead of fake positions fakePositions.insert(fakePositions.end(), seg.fakePositions.begin(), seg.fakePositions.end()); for (size_t i = 0; i < seg.clusterSet.size(); ++i) { if (seg.clusterSet[i]) clusterSet.set(i); } } + + if (segments.size() > 1) { + direction = positions[1] - positions[0]; + } else { + direction = segments[0].momentum; + } + + std::vector<Acts::Vector2> pointsZX {}; + for (const Acts::Vector3 &pos : fakePositions) { + pointsZX.push_back({pos.z(), pos.x()}); + } + linearFit(pointsZX); + + if (segments.size() > 1) { + fakeFit(); + } else { + momentum = 9999999.; + charge = 1; + } + getChi2(); size = clusters.size(); } @@ -253,9 +282,8 @@ CircleFitTrackSeedTool::Seed::Seed(const std::vector<Segment> &segments) : void CircleFitTrackSeedTool::Seed::fit() { CircleFit::CircleData circleData(positions); CircleFit::Circle circle = CircleFit::circleFit(circleData); - momentum = circle.r * 0.001 * 0.3 * 0.55; - charge = circle.cy > 0 ? 1 : -1; - // std::cout << "circle: \nnp.array([" << circle.cx << ", " << circle.cy << ", " << circle.r << "])" << std::endl; + momentum = circle.r > 0 ? circle.r * 0.001 * 0.3 * 0.55 : 9999999.; + charge = circle.cy < 0 ? -1 : 1; } void CircleFitTrackSeedTool::Seed::fakeFit(double B) { @@ -268,8 +296,8 @@ void CircleFitTrackSeedTool::Seed::fakeFit(double B) { charge = circle.cy > 0 ? 1 : -1; } -void CircleFitTrackSeedTool::Seed::linearFit() { - auto [origin, dir] = LinearFit::linearFit(positions); +void CircleFitTrackSeedTool::Seed::linearFit(const std::vector<Acts::Vector2> &points) { + auto [origin, dir] = LinearFit::linearFit(points); c1 = dir[1]/dir[0]; c0 = origin[1] - origin[0] * c1; } @@ -286,9 +314,6 @@ double CircleFitTrackSeedTool::Seed::getX(double z) { } void CircleFitTrackSeedTool::Seed::getChi2() { - linearFit(); - fakeFit(); - chi2 = 0; for (const Acts::Vector3 &pos : fakePositions) { dy = pos.y() - getY(pos.z()); @@ -302,9 +327,7 @@ void CircleFitTrackSeedTool::Seed::getChi2() { } Acts::CurvilinearTrackParameters CircleFitTrackSeedTool::Seed::get_params(double origin, Acts::BoundSymMatrix cov) const { - Acts::Vector3 dir = positions[1] - positions[0]; - dir = dir.normalized(); - Acts::Vector3 pos = positions[0] - (positions[0].z() - origin)/dir.z() * dir; + Acts::Vector3 pos = positions[0] - (positions[0].z() - origin)/direction.z() * direction; Acts::Vector4 pos4 {pos.x(), pos.y(), pos.z(), 0}; - return Acts::CurvilinearTrackParameters(pos4, dir, momentum, charge, cov); + return Acts::CurvilinearTrackParameters(pos4, direction.normalized(), momentum, charge, cov); } diff --git a/Tracking/Acts/FaserActsKalmanFilter/src/KalmanFitterTool.cxx b/Tracking/Acts/FaserActsKalmanFilter/src/KalmanFitterTool.cxx index de5f16132912d82827b21ece5c8c13297cac6252..6771d30cd26fa062b5ff7a6e080f00ed9dd69dc3 100644 --- a/Tracking/Acts/FaserActsKalmanFilter/src/KalmanFitterTool.cxx +++ b/Tracking/Acts/FaserActsKalmanFilter/src/KalmanFitterTool.cxx @@ -286,7 +286,7 @@ KalmanFitterTool::makeTrack(Acts::GeometryContext& geoCtx, TrackFitterResult& fi state.smoothed(), state.smoothedCovariance()); actsSmoothedParam.push_back(std::make_unique<const Acts::BoundTrackParameters>(Acts::BoundTrackParameters(actsParam))); // const auto& psurface=actsParam.referenceSurface(); - Acts::Vector2 local(actsParam.parameters()[Acts::eBoundLoc0], actsParam.parameters()[Acts::eBoundLoc1]); + // Acts::Vector2 local(actsParam.parameters()[Acts::eBoundLoc0], actsParam.parameters()[Acts::eBoundLoc1]); // const Acts::Vector3 dir = Acts::makeDirectionUnitFromPhiTheta(actsParam.parameters()[Acts::eBoundPhi], actsParam.parameters()[Acts::eBoundTheta]); // auto pos=actsParam.position(tgContext); parm = ConvertActsTrackParameterToATLAS(actsParam, geoCtx); @@ -310,7 +310,7 @@ KalmanFitterTool::makeTrack(Acts::GeometryContext& geoCtx, TrackFitterResult& fi double nDoF = state.calibratedSize(); const Trk::FitQualityOnSurface *quality = new Trk::FitQualityOnSurface(state.chi2(), nDoF); const Trk::TrackStateOnSurface *perState = new Trk::TrackStateOnSurface(measState, parm, quality, nullptr, typePattern); - // If a state was succesfully created add it to the trajectory + // If a state was successfully created add it to the trajectory if (perState) { finalTrajectory->insert(finalTrajectory->begin(), perState); } diff --git a/Tracking/Acts/FaserActsKalmanFilter/src/RootTrajectoryStatesWriterTool.cxx b/Tracking/Acts/FaserActsKalmanFilter/src/RootTrajectoryStatesWriterTool.cxx index 3964e0f0b98796ff840ee5597da6a8a75d0901c9..0444960bf295ddc05743ac8c2c047bec7e477fb5 100644 --- a/Tracking/Acts/FaserActsKalmanFilter/src/RootTrajectoryStatesWriterTool.cxx +++ b/Tracking/Acts/FaserActsKalmanFilter/src/RootTrajectoryStatesWriterTool.cxx @@ -243,7 +243,7 @@ StatusCode RootTrajectoryStatesWriterTool::write(const Acts::GeometryContext& gc simData = std::make_shared<TrackerSimDataCollection>(*simDataHandle); SG::ReadHandle<FaserSiHitCollection> siHitCollection {m_faserSiHitKey, ctx}; - ATH_CHECK(siHitCollection.isValid()); + ATH_CHECK(siHitCollection.isValid()); for (const FaserSiHit& hit : *siHitCollection) { int barcode = hit.trackNumber(); Identifier id = m_idHelper->wafer_id(hit.getStation(), hit.getPlane(), hit.getRow(), hit.getModule(), hit.getSensor()); diff --git a/Tracking/Acts/FaserActsKalmanFilter/src/RootTrajectorySummaryWriterTool.cxx b/Tracking/Acts/FaserActsKalmanFilter/src/RootTrajectorySummaryWriterTool.cxx index b4b97c1586048f027a7d7e353e4da88b8df2cfb9..e6afe3a77cf95528ac17d94b89c019048a7876a8 100644 --- a/Tracking/Acts/FaserActsKalmanFilter/src/RootTrajectorySummaryWriterTool.cxx +++ b/Tracking/Acts/FaserActsKalmanFilter/src/RootTrajectorySummaryWriterTool.cxx @@ -118,9 +118,11 @@ StatusCode RootTrajectorySummaryWriterTool::initialize() { StatusCode RootTrajectorySummaryWriterTool::finalize() { - m_outputFile->cd(); - m_outputTree->Write(); - m_outputFile->Close(); + if (!m_noDiagnostics) { + m_outputFile->cd(); + m_outputTree->Write(); + m_outputFile->Close(); + } return StatusCode::SUCCESS; } diff --git a/Tracking/Acts/FaserActsKalmanFilter/test/CKF2.py b/Tracking/Acts/FaserActsKalmanFilter/test/CKF2.py index 6e65c2e148d1a890b8f363c47fa984bd1764b3a7..52605c20f3a278203e721d59674b70e0d0dac63c 100644 --- a/Tracking/Acts/FaserActsKalmanFilter/test/CKF2.py +++ b/Tracking/Acts/FaserActsKalmanFilter/test/CKF2.py @@ -7,7 +7,7 @@ from AthenaCommon.Configurable import Configurable from CalypsoConfiguration.AllConfigFlags import ConfigFlags from CalypsoConfiguration.MainServicesConfig import MainServicesCfg from AthenaPoolCnvSvc.PoolReadConfig import PoolReadCfg -# from AthenaPoolCnvSvc.PoolWriteConfig import PoolWriteCfg +from AthenaPoolCnvSvc.PoolWriteConfig import PoolWriteCfg from TrackerPrepRawDataFormation.TrackerPrepRawDataFormationConfig import FaserSCT_ClusterizationCfg from TrackerSpacePointFormation.TrackerSpacePointFormationConfig import TrackerSpacePointFinderCfg from TrackerSegmentFit.TrackerSegmentFitConfig import SegmentFitAlgCfg @@ -17,24 +17,27 @@ from FaserActsKalmanFilter.CKF2Config import CKF2Cfg log.setLevel(DEBUG) Configurable.configurableRun3Behavior = True -ConfigFlags.Input.Files = ['my.RDO.pool.root'] +ConfigFlags.Input.Files = ['/home/tboeckh/Documents/data/moc/rdo/101302/FaserMC-MDC_PG_mumi_logE-101302-00000-RDO.root'] +# ConfigFlags.Input.Files = ['my.RDO.pool.root'] +ConfigFlags.addFlag("Output.xAODFileName", f"CKF.xAOD.root") ConfigFlags.Output.ESDFileName = "CKF.ESD.pool.root" -ConfigFlags.IOVDb.GlobalTag = "OFLCOND-FASER-01" -ConfigFlags.GeoModel.FaserVersion = "FASER-01" +ConfigFlags.IOVDb.GlobalTag = "OFLCOND-FASER-02" +ConfigFlags.GeoModel.FaserVersion = "FASERNU-03" ConfigFlags.GeoModel.Align.Dynamic = False ConfigFlags.Beam.NumberOfCollisions = 0. ConfigFlags.TrackingGeometry.MaterialSource = "Input" +ConfigFlags.Input.isMC = True ConfigFlags.lock() acc = MainServicesCfg(ConfigFlags) acc.merge(PoolReadCfg(ConfigFlags)) -# acc.merge(PoolWriteCfg(ConfigFlags)) +acc.merge(PoolWriteCfg(ConfigFlags)) acc.merge(FaserSCT_ClusterizationCfg(ConfigFlags)) acc.merge(TrackerSpacePointFinderCfg(ConfigFlags)) acc.merge(SegmentFitAlgCfg(ConfigFlags, SharedHitFraction=0.61, MinClustersPerFit=5, TanThetaXZCut=0.083)) acc.merge(GhostBustersCfg(ConfigFlags)) -acc.merge(CKF2Cfg(ConfigFlags, noDiagnostics=False)) +acc.merge(CKF2Cfg(ConfigFlags, noDiagnostics=True)) acc.getEventAlgo("CKF2").OutputLevel = DEBUG # logging.getLogger('forcomps').setLevel(VERBOSE) @@ -45,5 +48,17 @@ acc.getEventAlgo("CKF2").OutputLevel = DEBUG # acc.printConfig(withDetails=True) # ConfigFlags.dump() +from OutputStreamAthenaPool.OutputStreamConfig import OutputStreamCfg +itemList = [ + "xAOD::EventInfo#*", + "xAOD::EventAuxInfo#*", + "xAOD::FaserTriggerData#*", + "xAOD::FaserTriggerDataAux#*", + "FaserSCT_RDO_Container#*", + "Tracker::FaserSCT_ClusterContainer#*", + "TrackCollection#*", +] +acc.merge(OutputStreamCfg(ConfigFlags, "xAOD", itemList)) + sc = acc.run(maxEvents=-1) sys.exit(not sc.isSuccess()) diff --git a/Tracking/Acts/FaserActsKalmanFilter/test/TI12_CKF.py b/Tracking/Acts/FaserActsKalmanFilter/test/TI12_CKF.py deleted file mode 100644 index 98acb56f00ace95d54277726768250e9543b06ab..0000000000000000000000000000000000000000 --- a/Tracking/Acts/FaserActsKalmanFilter/test/TI12_CKF.py +++ /dev/null @@ -1,56 +0,0 @@ -#!/usr/bin/env python - -import sys -from AthenaCommon.Logging import log, logging -from AthenaCommon.Constants import DEBUG, VERBOSE, INFO -from AthenaCommon.Configurable import Configurable -from CalypsoConfiguration.AllConfigFlags import ConfigFlags -from CalypsoConfiguration.MainServicesConfig import MainServicesCfg -from FaserByteStreamCnvSvc.FaserByteStreamCnvSvcConfig import FaserByteStreamCnvSvcCfg -from TrackerPrepRawDataFormation.TrackerPrepRawDataFormationConfig import FaserSCT_ClusterizationCfg -from TrackerSegmentFit.TrackerSegmentFitConfig import SegmentFitAlgCfg -from FaserActsKalmanFilter.CombinatorialKalmanFilterConfig import CombinatorialKalmanFilterCfg - -log.setLevel(DEBUG) -Configurable.configurableRun3Behavior = True - -ConfigFlags.Input.Files = ['/home/tboeckh/tmp/Faser-Physics-006470-00093.raw_middleStation.SPs'] -ConfigFlags.Output.ESDFileName = "CKF.ESD.pool.root" -ConfigFlags.IOVDb.GlobalTag = "OFLCOND-FASER-02" -ConfigFlags.GeoModel.FaserVersion = "FASER-02" -ConfigFlags.IOVDb.DatabaseInstance = "OFLP200" -ConfigFlags.Input.ProjectName = "data21" -ConfigFlags.Input.isMC = False -ConfigFlags.Common.isOnline = False -ConfigFlags.GeoModel.Align.Dynamic = False -ConfigFlags.Beam.NumberOfCollisions = 0. -ConfigFlags.Detector.GeometryFaserSCT = True -ConfigFlags.lock() - -acc = MainServicesCfg(ConfigFlags) -acc.merge(FaserByteStreamCnvSvcCfg(ConfigFlags)) -acc.merge(FaserSCT_ClusterizationCfg(ConfigFlags, name="LevelClustering", DataObjectName="SCT_LEVELMODE_RDOs", ClusterToolTimingPattern="X1X")) -acc.merge(SegmentFitAlgCfg(ConfigFlags, name=f"LevelFit", MaxClusters=44)) -# acc.merge(FaserSCT_ClusterizationCfg(ConfigFlags, DataObjectName="SCT_RDOs", ClusterToolTimingPattern="XXX")) -# acc.merge(SegmentFitAlgCfg(ConfigFlags)) -acc.merge(CombinatorialKalmanFilterCfg(ConfigFlags)) -acc.getEventAlgo("CombinatorialKalmanFilterAlg").OutputLevel = VERBOSE - -replicaSvc = acc.getService("DBReplicaSvc") -replicaSvc.COOLSQLiteVetoPattern = "" -replicaSvc.UseCOOLSQLite = True -replicaSvc.UseCOOLFrontier = False -replicaSvc.UseGeomSQLite = True - -logging.getLogger('forcomps').setLevel(VERBOSE) -acc.foreach_component("*").OutputLevel = VERBOSE -acc.foreach_component("*ClassID*").OutputLevel = INFO -acc.getCondAlgo("FaserSCT_AlignCondAlg").OutputLevel = VERBOSE -acc.getCondAlgo("FaserSCT_DetectorElementCondAlg").OutputLevel = VERBOSE -acc.getService("StoreGateSvc").Dump = True -acc.getService("ConditionStore").Dump = True -acc.printConfig(withDetails=True) -ConfigFlags.dump() - -sc = acc.run(maxEvents=-1) -sys.exit(not sc.isSuccess())