File indexing completed on 2024-08-06 04:43:02
0001 #include "FWCore/TestProcessor/interface/TestProcessor.h"
0002 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0003 #include "RecoVertex/BeamSpotProducer/plugins/BeamSpotCompatibilityChecker.cc"
0004
0005 #define CATCH_CONFIG_MAIN
0006 #include "catch.hpp"
0007
0008
0009 TEST_CASE("Significance Calculation", "[Significance]") {
0010 double a = 10.0;
0011 double b = 12.0;
0012 double errA = 1.0;
0013 double errB = 1.5;
0014
0015 Significance sig(a, b, errA, errB);
0016 float significance = sig.getSig(false);
0017
0018
0019 REQUIRE(significance == Approx(1.1094).epsilon(10e-6));
0020 }
0021
0022
0023 TEST_CASE("BeamSpot Compatibility Checker", "[compareBS]") {
0024 reco::BeamSpot::Point pointA(1.0, 2.0, 3.0);
0025 reco::BeamSpot::Point pointB(1.1, 2.1, 3.1);
0026
0027 reco::BeamSpot::CovarianceMatrix matrixA;
0028 reco::BeamSpot::CovarianceMatrix matrixB;
0029
0030
0031 for (int i = 0; i < reco::BeamSpot::dimension; ++i) {
0032 for (int j = 0; j < reco::BeamSpot::dimension; ++j) {
0033 matrixA(i, j) = 0.01 * (i + 1) * (j + 1);
0034 matrixB(i, j) = 0.02 * (i + 1) * (j + 1);
0035 }
0036 }
0037
0038 reco::BeamSpot beamSpotA(pointA, 4.0, 0.01, 0.01, 0.1, matrixA);
0039 reco::BeamSpot beamSpotB(pointB, 4.2, 0.01, 0.01, 0.12, matrixB);
0040
0041
0042 edm::ParameterSet pset;
0043 pset.addParameter<double>("warningThr", 1.0);
0044 pset.addParameter<double>("errorThr", 2.0);
0045 pset.addParameter<edm::InputTag>("bsFromEvent", edm::InputTag(""));
0046 pset.addParameter<bool>("dbFromEvent", false);
0047 pset.addParameter<edm::InputTag>("bsFromDB", edm::InputTag(""));
0048
0049 BeamSpotCompatibilityChecker checker(pset);
0050 auto significances = checker.compareBS(beamSpotA, beamSpotB);
0051
0052
0053 for (size_t i = 0; i < significances.size(); ++i) {
0054 std::cout << "Significance[" << i << "]: " << significances[i] << std::endl;
0055 }
0056
0057 REQUIRE(significances[0] == Approx(0.57735).epsilon(10e-6));
0058 REQUIRE(significances[1] == Approx(0.288675).epsilon(10e-6));
0059 REQUIRE(significances[2] == Approx(0.19245).epsilon(10e-6));
0060 REQUIRE(significances[3] == Approx(0.0164957).epsilon(10e-6));
0061 REQUIRE(significances[4] == Approx(0.0164957).epsilon(10e-6));
0062 REQUIRE(significances[5] == Approx(0.288675).epsilon(10e-6));
0063 }