tests/algorithms_test/reco_FarForwardNeutronReconstruction.cc

0001 // SPDX-License-Identifier: LGPL-3.0-or-later
0002 // Copyright (C) 2024, Sebouh Paul
0003
0004 #include <Evaluator/DD4hepUnits.h>                 // for MeV, mm, keV, ns
0005 #include <catch2/catch_test_macros.hpp>            // for AssertionHandler, operator""_catch_sr, StringRef, REQUIRE, operator<, operator==, operator>, TEST_CASE
0006 #include <edm4eic/ClusterCollection.h>
0007 #include <edm4eic/ReconstructedParticleCollection.h>
0008 #include <edm4hep/Vector3f.h>                      // for Vector3f
0009 #include <spdlog/common.h>                         // for level_enum
0010 #include <spdlog/logger.h>                         // for logger
0011 #include <spdlog/spdlog.h>                         // for default_logger
0012 #include <stddef.h>                                // for size_t
0013 #include <stdlib.h>
0014 #include <array>                                   // for array
0015 #include <cmath>                                   // for sqrt, abs
0016 #include <iostream>
0017 #include <memory>                                  // for allocator, unique_ptr, make_unique, shared_ptr, __shared_ptr_access
0018 #include <vector>
0019
0020 #include "algorithms/reco/FarForwardNeutronReconstruction.h"
0021 #include "algorithms/reco/FarForwardNeutronReconstructionConfig.h"
0022
0023 using eicrecon::FarForwardNeutronReconstruction;
0024 using eicrecon::FarForwardNeutronReconstructionConfig;
0025
0026 TEST_CASE( "the cluster merging algorithm runs", "[FarForwardNeutronReconstruction]" ) {
0027   FarForwardNeutronReconstruction algo("FarForwardNeutronReconstruction");
0028
0029   std::shared_ptr<spdlog::logger> logger = spdlog::default_logger()->clone("FarForwardNeutronReconstruction");
0030   logger->set_level(spdlog::level::trace);
0031
0032   FarForwardNeutronReconstructionConfig cfg;
0033   std::vector<double> corr_parameters={-0.0756, -1.91,  2.30};
0034   cfg.scale_corr_coeff_hcal=corr_parameters;
0035   cfg.scale_corr_coeff_ecal=corr_parameters;
0036   algo.applyConfig(cfg);
0037   algo.init();
0038
0039   edm4eic::ClusterCollection clust_coll_hcal;
0040   std::array<float,3> x={30*dd4hep::mm,90*dd4hep::mm,0};
0041   std::array<float,3> y={-30*dd4hep::mm,0*dd4hep::mm, -90*dd4hep::mm};
0042   std::array<float,3> z={30*dd4hep::m,30*dd4hep::m, 30*dd4hep::m};
0043   std::array<double,3> E={80*dd4hep::GeV,5*dd4hep::GeV,5*dd4hep::GeV};
0044   float sumEnergies=0;
0045   for(size_t i=0; i<3; i++){
0046     auto cluster=clust_coll_hcal.create();
0047     cluster.setEnergy(E[i]);
0048     cluster.setPosition({x[i], y[i], z[i]});
0049
0050   }
0051
0052   edm4eic::ClusterCollection clust_coll_ecal;
0053   auto ecal_cluster=clust_coll_ecal.create();
0054   ecal_cluster.setEnergy(2);
0055   ecal_cluster.setPosition({0, 0, 25*dd4hep::m});
0056
0057   auto neutroncand_coll = std::make_unique<edm4eic::ReconstructedParticleCollection>();
0058   algo.process({&clust_coll_hcal, &clust_coll_ecal}, {neutroncand_coll.get()});
0059
0060   REQUIRE( (*neutroncand_coll).size() == 1);
0061
0062   double corr=algo.calc_corr(92, corr_parameters);
0063   double tol=0.001;
0064   double E_expected=92*dd4hep::GeV*1/(1+corr);
0065   double Px_expected=0.09199*1/(1+corr);
0066   double Py_expected=-0.09199*1/(1+corr);
0067   double Pz_expected=91.99*1/(1+corr);
0068   //check that the correct energy and momenta are being obtained
0069   std::cout << "E, px, py, pz = " << (*neutroncand_coll)[0].getEnergy() <<"  " << (*neutroncand_coll)[0].getMomentum().x << "  " << (*neutroncand_coll)[0].getMomentum().y << "  " << (*neutroncand_coll)[0].getMomentum().z << std::endl;
0070   REQUIRE( abs((*neutroncand_coll)[0].getEnergy()-E_expected)/E_expected<tol);
0071   REQUIRE( abs((*neutroncand_coll)[0].getMomentum().x-Px_expected)/Px_expected<tol);
0072   REQUIRE( abs((*neutroncand_coll)[0].getMomentum().y-Py_expected)/Py_expected<tol);
0073   REQUIRE( abs((*neutroncand_coll)[0].getMomentum().z-Pz_expected)/Pz_expected<tol);
0074
0075 }