modules.cc
Go to the documentation of this file.
1 /*
2  * MoMEMta: a modular implementation of the Matrix Element Method
3  * Copyright (C) 2016 Universite catholique de Louvain (UCL), Belgium
4  *
5  * This program is free software: you can redistribute it and/or modify
6  * it under the terms of the GNU General Public License as published by
7  * the Free Software Foundation, either version 3 of the License, or
8  * (at your option) any later version.
9  *
10  * This program is distributed in the hope that it will be useful,
11  * but WITHOUT ANY WARRANTY; without even the implied warranty of
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13  * GNU General Public License for more details.
14  *
15  * You should have received a copy of the GNU General Public License
16  * along with this program. If not, see <http://www.gnu.org/licenses/>.
17  */
18 
26 #include <catch.hpp>
27 
28 #include <momemta/Configuration.h>
29 #include <momemta/ModuleFactory.h>
30 #include <momemta/Module.h>
31 #include <momemta/ParameterSet.h>
32 #include <momemta/Pool.h>
33 #include <momemta/Solution.h>
34 #include <momemta/Types.h>
35 #include <momemta/Math.h>
36 
37 #define N_PS_POINTS 5
38 
39 // A mock of ParameterSet to change visibility of the constructor
41 public:
42  ParameterSetMock(const std::string& name):
43  ParameterSet(name, name) {
44  // Empty
45  }
46 
47  void createMock(const std::string& name, const momemta::any& value) {
48  create(name, value);
49  }
50 };
51 
52 std::shared_ptr<std::vector<double>> addPhaseSpacePoints(std::shared_ptr<Pool> pool) {
53  auto ps_points = pool->put<std::vector<double>>({"cuba", "ps_points"});
54  auto weight = pool->put<double>({"cuba", "ps_weight"});
55 
56  double step = 1. / (N_PS_POINTS - 1);
57 
58  ps_points->resize(N_PS_POINTS);
59  for (size_t i = 0; i < N_PS_POINTS; i++) {
60  double point = i * step;
61  ps_points->operator[](i) = point;
62  }
63 
64  *weight = 1;
65 
66  return ps_points;
67 }
68 
69 std::shared_ptr<std::vector<LorentzVector>> addInputParticles(std::shared_ptr<Pool> pool) {
70  auto inputs = pool->put<std::vector<LorentzVector>>({"input", "particles"});
71 
72  // Some random massive 4-vectors to be used by the tests
73  inputs->push_back( { 16.171895980835, -13.7919054031372, -3.42997527122497, 21.5293197631836 });
74  inputs->push_back( { 71.3899612426758, 96.0094833374023, -77.2513122558594, 142.492813110352 });
75  inputs->push_back( { -18.9018573760986, 10.0896110534668, -0.602926552295686, 21.4346446990967 });
76  inputs->push_back( { -55.7908325195313, -111.59294128418, -122.144721984863, 174.66259765625 });
77 
78  // Some random massless 4-vectors to be used by the tests
79  inputs->push_back( { 26.2347,76.1854,32.6172,86.9273 } );
80  inputs->push_back( { -27.3647,30.2536,38.5939,56.1569 } );
81  inputs->push_back( { -7.0817,-63.9966,-34.8497,73.2135 } );
82  inputs->push_back( { 80.756,3.12662,33.0355,87.3078 } );
83 
84  return inputs;
85 }
86 
87 TEST_CASE("Modules", "[modules]") {
88  std::shared_ptr<Pool> pool(new Pool());
89  std::shared_ptr<ParameterSetMock> parameters;
90 
91  auto createModule = [&pool, &parameters](const std::string& type) {
93  module.name = type;
94  module.type = type;
95  module.parameters.reset(parameters->clone());
96 
97  auto result = momemta::ModuleFactory::get().create(type, pool, *parameters);
98  REQUIRE(result.get());
99 
100  return result;
101  };
102 
103  // Register phase space points, mocking what cuba would do
104  auto ps_points = addPhaseSpacePoints(pool);
105  // Put a few random input particles into the pool
106  auto input_particles = addInputParticles(pool);
107 
108  SECTION("BreitWignerGenerator") {
109 
110  parameters.reset(new ParameterSetMock("BreitWignerGenerator"));
111 
112  double mass = 173.;
113  double width = 5.;
114 
115  parameters->set("mass", mass);
116  parameters->set("width", width);
117  parameters->set("ps_point", InputTag("cuba", "ps_points", 2)); // 0.5 by default
118 
119  auto module = createModule("BreitWignerGenerator");
120 
121  auto s = pool->get<double>({"BreitWignerGenerator", "s"});
122  auto jacobian = pool->get<double>({"BreitWignerGenerator", "jacobian"});
123 
124  REQUIRE(module->work() == Module::Status::OK);
125 
126  REQUIRE(*s == Approx(29941.5));
127  REQUIRE(*jacobian == Approx(2693.05));
128 
129  ps_points->operator[](2) = 0;
130 
131  REQUIRE(module->work() == Module::Status::OK);
132 
133  REQUIRE(*s == Approx(0).margin(std::numeric_limits<float>::epsilon()));
134  }
135 
136  SECTION("UniformGenerator") {
137 
138  parameters.reset(new ParameterSetMock("UniformGenerator"));
139 
140  double min = 100;
141  double max = 250;
142  double expected_jacobian = max - min;
143 
144  parameters->set("min", min);
145  parameters->set("max", max);
146  parameters->set("ps_point", InputTag("cuba", "ps_points", 0));
147 
148  auto module = createModule("UniformGenerator");
149 
150  auto output = pool->get<double>({"UniformGenerator", "output"});
151  auto jacobian = pool->get<double>({"UniformGenerator", "jacobian"});
152 
153  ps_points->operator[](0) = 0;
154 
155  REQUIRE(module->work() == Module::Status::OK);
156 
157  REQUIRE(*output == Approx(min));
158  REQUIRE(*jacobian == Approx(expected_jacobian));
159 
160  ps_points->operator[](0) = 1;
161 
162  REQUIRE(module->work() == Module::Status::OK);
163 
164  REQUIRE(*output == Approx(max));
165  REQUIRE(*jacobian == Approx(expected_jacobian));
166 
167  ps_points->operator[](0) = 0.5;
168 
169  REQUIRE(module->work() == Module::Status::OK);
170 
171  REQUIRE(*output == Approx((max + min) / 2.));
172  REQUIRE(*jacobian == Approx(expected_jacobian));
173  }
174 
175  SECTION("BlockA") {
176 
177  parameters.reset(new ParameterSetMock("BlockA"));
178 
179  double sqrt_s = 13000;
180 
181  parameters->set("energy", sqrt_s);
182 
183  parameters->set("p1", InputTag("input", "particles", 0));
184  parameters->set("p2", InputTag("input", "particles", 1));
185  parameters->set("branches", std::vector<InputTag>( { InputTag("input", "particles", 2) } ));
186 
187  Value<SolutionCollection> solutions = pool->get<SolutionCollection>({"BlockA", "solutions"});
188 
189  auto module = createModule("BlockA");
190 
191  REQUIRE(module->work() == Module::Status::OK);
192  REQUIRE(solutions->size() == 1);
193 
194  for (const auto& solution: *solutions) {
195  REQUIRE(solution.valid == true);
196 
197  LorentzVector test_pT = solution.values.at(0) + solution.values.at(1) + input_particles->at(2);
198  REQUIRE(test_pT.Pt() == Approx(0).margin(std::numeric_limits<float>::epsilon()));
199 
200  REQUIRE(solution.values.at(0).M() == Approx(input_particles->at(0).M()));
201  REQUIRE(solution.values.at(0).Phi() == Approx(input_particles->at(0).Phi()));
202  REQUIRE(solution.values.at(0).Theta() == Approx(input_particles->at(0).Theta()));
203 
204  REQUIRE(solution.values.at(1).M() == Approx(input_particles->at(1).M()));
205  REQUIRE(solution.values.at(1).Phi() == Approx(input_particles->at(1).Phi()));
206  REQUIRE(solution.values.at(1).Theta() == Approx(input_particles->at(1).Theta()));
207  }
208  }
209 
210  SECTION("BlockB") {
211 
212  parameters.reset(new ParameterSetMock("BlockB"));
213 
214  double sqrt_s = 13000;
215  bool pT_is_met = false;
216 
217  parameters->set("energy", sqrt_s);
218  parameters->set("pT_is_met", pT_is_met);
219 
220  parameters->set("s12", InputTag("mockS", "s12"));
221  double s_12 = SQ(200);
222  auto s12 = pool->put<double>({"mockS", "s12"});
223  *s12 = s_12;
224 
225  parameters->set("p2", InputTag("input", "particles", 0));
226 
227  Value<SolutionCollection> solutions = pool->get<SolutionCollection>({"BlockB", "solutions"});
228 
229  auto module = createModule("BlockB");
230 
231  REQUIRE(module->work() == Module::Status::OK);
232  REQUIRE(solutions->size() == 2);
233 
234  for (const auto& solution: *solutions) {
235  REQUIRE(solution.valid == true);
236 
237  LorentzVector test_p12 = input_particles->at(0) + solution.values.at(0);
238 
239  REQUIRE(test_p12.M2() == Approx(s_12));
240  REQUIRE(test_p12.Pt() == Approx(0));
241  REQUIRE(solution.values.at(0).M() / solution.values.at(0).E() == Approx(0).margin(std::numeric_limits<float>::epsilon()));
242  }
243  }
244 
245  SECTION("BlockC") {
246 
247  parameters.reset(new ParameterSetMock("BlockC"));
248 
249  double sqrt_s = 13000;
250  bool pT_is_met = false;
251  double m1 = 1.;
252 
253  parameters->set("energy", sqrt_s);
254  parameters->set("pT_is_met", pT_is_met);
255  parameters->set("m1", m1);
256 
257  parameters->set("s12", InputTag("mockS", "s12"));
258  parameters->set("s123", InputTag("mockS", "s123"));
259 
260  double s_12 = SQ(80);
261  double s_123 = SQ(173);
262 
263  auto s12 = pool->put<double>({"mockS", "s12"});
264  auto s123 = pool->put<double>({"mockS", "s123"});
265 
266  *s12 = s_12;
267  *s123 = s_123;
268 
269  parameters->set("p2", InputTag("input", "particles", 0));
270  parameters->set("p3", InputTag("input", "particles", 5));
271  parameters->set("branches", std::vector<InputTag>({InputTag("input", "particles", 2)}));
272 
273  Value<SolutionCollection> solutions = pool->get<SolutionCollection>({"BlockC", "solutions"});
274 
275  auto module = createModule("BlockC");
276 
277  REQUIRE(module->work() == Module::Status::OK);
278  REQUIRE(solutions->size() >= 2);
279 
280  for (const auto& solution : *solutions) {
281  REQUIRE(solution.valid == true);
282 
283  LorentzVector test_p12 = input_particles->at(0) + solution.values.at(0);
284  LorentzVector test_p123 = solution.values.at(1) + test_p12;
285 
286  REQUIRE(test_p12.M2() == Approx(s_12));
287  REQUIRE(test_p123.M2() == Approx(s_123));
288 
289  LorentzVector test_pT = test_p123 + input_particles->at(2);
290  REQUIRE(test_pT.Pt() == Approx(0).margin(std::numeric_limits<float>::epsilon()));
291 
292  REQUIRE(solution.values.at(0).M() == Approx(m1));
293  REQUIRE(solution.values.at(1).M() / solution.values.at(1).E() == Approx(0).margin(std::numeric_limits<float>::epsilon()));
294  }
295  }
296 
297  SECTION("BlockD") {
298 
299  parameters.reset(new ParameterSetMock("BlockD"));
300 
301  double sqrt_s = 13000;
302  bool pT_is_met = false;
303 
304  parameters->set("energy", sqrt_s);
305  parameters->set("pT_is_met", pT_is_met);
306 
307  parameters->set("s13", InputTag("mockS", "s13"));
308  parameters->set("s134", InputTag("mockS", "s134"));
309  parameters->set("s25", InputTag("mockS", "s25"));
310  parameters->set("s256", InputTag("mockS", "s256"));
311 
312  double s_13_25 = SQ(80);
313  double s_134_256 = SQ(170);
314 
315  auto s13 = pool->put<double>({"mockS", "s13"});
316  auto s134 = pool->put<double>({"mockS", "s134"});
317  auto s25 = pool->put<double>({"mockS", "s25"});
318  auto s256 = pool->put<double>({"mockS", "s256"});
319 
320  *s13 = s_13_25;
321  *s134 = s_134_256;
322  *s25 = s_13_25;
323  *s256 = s_134_256;
324 
325  parameters->set("p3", InputTag("input", "particles", 0));
326  parameters->set("p4", InputTag("input", "particles", 1));
327  parameters->set("p5", InputTag("input", "particles", 2));
328  parameters->set("p6", InputTag("input", "particles", 3));
329 
330  Value<SolutionCollection> solutions = pool->get<SolutionCollection>({"BlockD", "solutions"});
331 
332  auto module = createModule("BlockD");
333 
334  REQUIRE(module->work() == Module::Status::OK);
335  REQUIRE(solutions->size() == 2);
336 
337  for (const auto& solution: *solutions) {
338  REQUIRE(solution.valid == true);
339 
340  LorentzVector test_p13 = input_particles->at(0) + solution.values.at(0);
341  LorentzVector test_p134 = input_particles->at(1) + test_p13;
342  LorentzVector test_p25 = input_particles->at(2) + solution.values.at(1);
343  LorentzVector test_p256 = input_particles->at(3) + test_p25;
344 
345  REQUIRE(test_p13.M2() == Approx(s_13_25));
346  REQUIRE(test_p134.M2() == Approx(s_134_256));
347  REQUIRE(test_p25.M2() == Approx(s_13_25));
348  REQUIRE(test_p256.M2() == Approx(s_134_256));
349 
350  LorentzVector test_pT = test_p134 + test_p256;
351  REQUIRE(test_pT.Pt() == Approx(0).margin(std::numeric_limits<float>::epsilon()));
352 
353  REQUIRE(solution.values.at(0).M() / solution.values.at(0).E() == Approx(0).margin(std::numeric_limits<float>::epsilon()));
354  REQUIRE(solution.values.at(1).M() / solution.values.at(1).E() == Approx(0).margin(std::numeric_limits<float>::epsilon()));
355  }
356  }
357 
358  SECTION("BlockE") {
359 
360  parameters.reset(new ParameterSetMock("BlockE"));
361 
362  double sqrt_s = 13000;
363  parameters->set("energy", sqrt_s);
364 
365  parameters->set("s13", InputTag("mockS", "s13"));
366  parameters->set("s24", InputTag("mockS", "s24"));
367  parameters->set("s_hat", InputTag("mockS", "s_hat"));
368  parameters->set("y_tot", InputTag("mockS", "y_tot"));
369 
370  double s_13 = SQ(100);
371  double s_24 = SQ(200);
372  double s_hat = SQ(500);
373  double rap = 1.5;
374 
375  *pool->put<double>({"mockS", "s13"}) = s_13;
376  *pool->put<double>({"mockS", "s24"}) = s_24;
377  *pool->put<double>({"mockS", "s_hat"}) = s_hat;
378  *pool->put<double>({"mockS", "y_tot"}) = rap;
379 
380  parameters->set("p3", InputTag("input", "particles", 0));
381  parameters->set("p4", InputTag("input", "particles", 2));
382 
383  Value<SolutionCollection> solutions = pool->get<SolutionCollection>({"BlockE", "solutions"});
384 
385  auto module = createModule("BlockE");
386 
387  REQUIRE(module->work() == Module::Status::OK);
388  REQUIRE(solutions->size() == 2);
389 
390  for (const auto& solution: *solutions) {
391  REQUIRE(solution.valid == true);
392 
393  LorentzVector test_p13 = input_particles->at(0) + solution.values.at(0);
394  LorentzVector test_p24 = input_particles->at(2) + solution.values.at(1);
395  LorentzVector test_pT = test_p13 + test_p24;
396 
397  REQUIRE(test_p13.M2() == Approx(s_13));
398  REQUIRE(test_p24.M2() == Approx(s_24));
399 
400  REQUIRE(test_pT.M2() == Approx(s_hat));
401  REQUIRE(test_pT.Rapidity() == Approx(rap));
402  REQUIRE(test_pT.Pt() == Approx(0));
403 
404  REQUIRE(solution.values.at(0).M() / solution.values.at(0).E() == Approx(0).margin(0.000001));
405  REQUIRE(solution.values.at(1).M() / solution.values.at(1).E() == Approx(0).margin(0.000001));
406  }
407  }
408 
409  SECTION("BlockF") {
410 
411  parameters.reset(new ParameterSetMock("BlockF"));
412 
413  double sqrt_s = 13000;
414  parameters->set("energy", sqrt_s);
415 
416  parameters->set("s13", InputTag("mockS", "s13"));
417  parameters->set("s24", InputTag("mockS", "s24"));
418  parameters->set("q1", InputTag("mockS", "q1"));
419  parameters->set("q2", InputTag("mockS", "q2"));
420 
421  double s_13 = SQ(100);
422  double s_24 = SQ(200);
423  double q1 = 0.172;
424  double q2 = 0.0086;
425 
426  *pool->put<double>({"mockS", "s13"}) = s_13;
427  *pool->put<double>({"mockS", "s24"}) = s_24;
428  *pool->put<double>({"mockS", "q1"}) = q1;
429  *pool->put<double>({"mockS", "q2"}) = q2;
430 
431  parameters->set("p3", InputTag("input", "particles", 0));
432  parameters->set("p4", InputTag("input", "particles", 2));
433 
434  Value<SolutionCollection> solutions = pool->get<SolutionCollection>({"BlockF", "solutions"});
435 
436  auto module = createModule("BlockF");
437 
438  REQUIRE(module->work() == Module::Status::OK);
439  REQUIRE(solutions->size() == 2);
440 
441  for (const auto& solution: *solutions) {
442  REQUIRE(solution.valid == true);
443 
444  LorentzVector test_p13 = input_particles->at(0) + solution.values.at(0);
445  LorentzVector test_p24 = input_particles->at(2) + solution.values.at(1);
446  LorentzVector test_pT = test_p13 + test_p24;
447 
448  REQUIRE(test_p13.M2() == Approx(s_13));
449  REQUIRE(test_p24.M2() == Approx(s_24));
450 
451  REQUIRE(test_pT.E() == Approx(0.5 * (q1 + q2) * sqrt_s));
452  REQUIRE(test_pT.Pz() == Approx(0.5 * (q1 - q2) * sqrt_s));
453  REQUIRE(test_pT.Pt() == Approx(0).margin(std::numeric_limits<float>::epsilon()));
454 
455  REQUIRE(solution.values.at(0).M() / solution.values.at(0).E() == Approx(0).margin(0.000001));
456  REQUIRE(solution.values.at(1).M() / solution.values.at(1).E() == Approx(0).margin(0.000001));
457  }
458  }
459 
460  SECTION("BlockG") {
461 
462  parameters.reset(new ParameterSetMock("BlockG"));
463 
464  double sqrt_s = 13000;
465  parameters->set("energy", sqrt_s);
466 
467  parameters->set("s12", InputTag("mockS", "s12"));
468  parameters->set("s34", InputTag("mockS", "s34"));
469 
470  double s_12 = SQ(650);
471  double s_34 = SQ(550);
472 
473  *pool->put<double>({"mockS", "s12"}) = s_12;
474  *pool->put<double>({"mockS", "s34"}) = s_34;
475 
476  parameters->set("p1", InputTag("input", "particles", 4));
477  parameters->set("p2", InputTag("input", "particles", 5));
478  parameters->set("p3", InputTag("input", "particles", 6));
479  parameters->set("p4", InputTag("input", "particles", 7));
480 
481  Value<SolutionCollection> solutions = pool->get<SolutionCollection>({"BlockG", "solutions"});
482 
483  auto module = createModule("BlockG");
484 
485  REQUIRE(module->work() == Module::Status::OK);
486  REQUIRE(solutions->size() == 1);
487 
488  for (const auto& solution: *solutions) {
489  REQUIRE(solution.valid == true);
490 
491  LorentzVector test_p12 = solution.values.at(0) + solution.values.at(1);
492  LorentzVector test_p34 = solution.values.at(2) + solution.values.at(3);
493  LorentzVector test_pT = test_p12 + test_p34;
494 
495  REQUIRE(test_p12.M2() == Approx(s_12));
496  REQUIRE(test_p34.M2() == Approx(s_34));
497 
498  REQUIRE(test_pT.Pt() == Approx(0).margin(std::numeric_limits<float>::epsilon()));
499 
500  REQUIRE(solution.values.at(0).M() / solution.values.at(0).E() == Approx(0).margin(std::numeric_limits<float>::epsilon()));
501  REQUIRE(solution.values.at(0).Phi() == Approx(input_particles->at(4).Phi()));
502  REQUIRE(solution.values.at(0).Theta() == Approx(input_particles->at(4).Theta()));
503 
504  REQUIRE(solution.values.at(1).M() / solution.values.at(1).E() == Approx(0).margin(std::numeric_limits<float>::epsilon()));
505  REQUIRE(solution.values.at(1).Phi() == Approx(input_particles->at(5).Phi()));
506  REQUIRE(solution.values.at(1).Theta() == Approx(input_particles->at(5).Theta()));
507 
508  REQUIRE(solution.values.at(2).M() / solution.values.at(2).E() == Approx(0).margin(std::numeric_limits<float>::epsilon()));
509  REQUIRE(solution.values.at(2).Phi() == Approx(input_particles->at(6).Phi()));
510  REQUIRE(solution.values.at(2).Theta() == Approx(input_particles->at(6).Theta()));
511 
512  REQUIRE(solution.values.at(3).M() / solution.values.at(3).E() == Approx(0).margin(std::numeric_limits<float>::epsilon()));
513  REQUIRE(solution.values.at(3).Phi() == Approx(input_particles->at(7).Phi()));
514  REQUIRE(solution.values.at(3).Theta() == Approx(input_particles->at(7).Theta()));
515 
516  }
517  }
518 
519  SECTION("SecondaryBlockA") {
520 
521  parameters.reset(new ParameterSetMock("SecondaryBlockA"));
522 
523  double sqrt_s = 13000;
524  parameters->set("energy", sqrt_s);
525 
526  parameters->set("s12", InputTag("mockS", "s12"));
527  parameters->set("s123", InputTag("mockS", "s123"));
528  parameters->set("s1234", InputTag("mockS", "s1234"));
529 
530  double s_12 = SQ(75);
531  double s_123 = SQ(150);
532  double s_1234 = SQ(350);
533 
534  *pool->put<double>({"mockS", "s12"}) = s_12;
535  *pool->put<double>({"mockS", "s123"}) = s_123;
536  *pool->put<double>({"mockS", "s1234"}) = s_1234;
537 
538  parameters->set("p2", InputTag("input", "particles", 0));
539  parameters->set("p3", InputTag("input", "particles", 2));
540  parameters->set("p4", InputTag("input", "particles", 1));
541 
542  Value<SolutionCollection> solutions = pool->get<SolutionCollection>({"SecondaryBlockA", "solutions"});
543 
544  auto module = createModule("SecondaryBlockA");
545 
546  REQUIRE(module->work() == Module::Status::OK);
547  REQUIRE(solutions->size() == 2);
548 
549  for (const auto& solution: *solutions) {
550  REQUIRE(solution.valid == true);
551 
552  LorentzVector test_p12 = input_particles->at(0) + solution.values.at(0);
553  LorentzVector test_p123 = input_particles->at(2) + test_p12;
554  LorentzVector test_p1234 = input_particles->at(1) + test_p123;
555 
556  REQUIRE(test_p12.M2() == Approx(s_12));
557  REQUIRE(test_p123.M2() == Approx(s_123));
558  REQUIRE(test_p1234.M2() == Approx(s_1234));
559 
560  REQUIRE(solution.values.at(0).M() / solution.values.at(0).E() == Approx(0).margin(0.000001));
561  }
562  }
563 
564  SECTION("SecondaryBlockB") {
565 
566  parameters.reset(new ParameterSetMock("SecondaryBlockB"));
567 
568  double sqrt_s = 13000;
569  parameters->set("energy", sqrt_s);
570 
571  parameters->set("s12", InputTag("mockS", "s12"));
572  parameters->set("s123", InputTag("mockS", "s123"));
573 
574  double s_12 = SQ(475.073);
575  double s_123 = SQ(516.951);
576 
577  *pool->put<double>({"mockS", "s12"}) = s_12;
578  *pool->put<double>({"mockS", "s123"}) = s_123;
579 
580  parameters->set("p1", InputTag("input", "particles", 0));
581  parameters->set("p2", InputTag("input", "particles", 1));
582  parameters->set("p3", InputTag("input", "particles", 2));
583 
584  Value<SolutionCollection> solutions = pool->get<SolutionCollection>({"SecondaryBlockB", "solutions"});
585 
586  auto module = createModule("SecondaryBlockB");
587 
588  REQUIRE(module->work() == Module::Status::OK);
589  REQUIRE(solutions->size() == 1);
590 
591  for (const auto& solution: *solutions) {
592  REQUIRE(solution.valid == true);
593 
594  LorentzVector test_p12 = input_particles->at(1) + solution.values.at(0);
595  LorentzVector test_p123 = input_particles->at(2) + test_p12;
596 
597  REQUIRE(test_p12.M2() == Approx(s_12));
598  REQUIRE(test_p123.M2() == Approx(s_123));
599 
600  REQUIRE(solution.values.at(0).M() == Approx(input_particles->at(0).M()));
601  REQUIRE(solution.values.at(0).Phi() == Approx(input_particles->at(0).Phi()));
602  }
603  }
604 
605  SECTION("SecondaryBlockCD") {
606 
607  parameters.reset(new ParameterSetMock("SecondaryBlockCD"));
608 
609  double sqrt_s = 13000;
610  parameters->set("energy", sqrt_s);
611 
612  parameters->set("s12", InputTag("mockS", "s12"));
613  double s_12 = SQ(50);
614  *pool->put<double>({"mockS", "s12"}) = s_12;
615 
616  parameters->set("p1", InputTag("input", "particles", 0));
617  parameters->set("p2", InputTag("input", "particles", 2));
618 
619  Value<SolutionCollection> solutions = pool->get<SolutionCollection>({"SecondaryBlockCD", "solutions"});
620 
621  auto module = createModule("SecondaryBlockCD");
622 
623  REQUIRE(module->work() == Module::Status::OK);
624  REQUIRE(solutions->size() == 1);
625 
626  for (const auto& solution: *solutions) {
627  REQUIRE(solution.valid == true);
628 
629  LorentzVector test_p12 = input_particles->at(2) + solution.values.at(0);
630 
631  REQUIRE(test_p12.M2() == Approx(s_12));
632 
633  REQUIRE(solution.values.at(0).Phi() == Approx(input_particles->at(0).Phi()));
634  REQUIRE(solution.values.at(0).Theta() == Approx(input_particles->at(0).Theta()));
635  }
636  }
637 
638  SECTION("SecondaryBlockE") {
639 
640  parameters.reset(new ParameterSetMock("SecondaryBlockE"));
641 
642  double sqrt_s = 13000;
643  parameters->set("energy", sqrt_s);
644 
645  parameters->set("s12", InputTag("mockS", "s12"));
646  parameters->set("s123", InputTag("mockS", "s123"));
647 
648  double s_12 = SQ(50);
649  double s_123 = SQ(500);
650 
651  *pool->put<double>({"mockS", "s12"}) = s_12;
652  *pool->put<double>({"mockS", "s123"}) = s_123;
653 
654  parameters->set("p1", InputTag("input", "particles", 0));
655  parameters->set("p2", InputTag("input", "particles", 5));
656  parameters->set("p3", InputTag("input", "particles", 3));
657 
658  Value<SolutionCollection> solutions = pool->get<SolutionCollection>({"SecondaryBlockE", "solutions"});
659 
660  auto module = createModule("SecondaryBlockE");
661 
662  REQUIRE(module->work() == Module::Status::OK);
663  REQUIRE(solutions->size() == 2);
664 
665  for (const auto& solution: *solutions) {
666  REQUIRE(solution.valid == true);
667 
668  LorentzVector test_p12 = solution.values.at(0) + solution.values.at(1);
669  LorentzVector test_p123 = input_particles->at(3) + test_p12;
670 
671  REQUIRE(test_p12.M2() == Approx(s_12));
672  REQUIRE(test_p123.M2() == Approx(s_123));
673 
674  REQUIRE(solution.values.at(0).M() == Approx(input_particles->at(0).M()));
675  REQUIRE(solution.values.at(0).Phi() == Approx(input_particles->at(0).Phi()));
676  REQUIRE(solution.values.at(0).Theta() == Approx(input_particles->at(0).Theta()));
677 
678  REQUIRE(solution.values.at(1).M() / solution.values.at(1).E() == Approx(0).margin(std::numeric_limits<float>::epsilon()));
679  REQUIRE(solution.values.at(1).Phi() == Approx(input_particles->at(5).Phi()));
680  REQUIRE(solution.values.at(1).Theta() == Approx(input_particles->at(5).Theta()));
681  }
682  }
683 }
Mathematical functions.
std::string name
Name of the module (user-defined from the configuration file)
Definition: Configuration.h:44
virtual void create(const std::string &name, const momemta::any &value)
Add a new element to the ParameterSet.
Definition: ParameterSet.cc:48
std::shared_ptr< ParameterSet > parameters
Module&#39;s parameters, as parsed from the configuration file.
Definition: Configuration.h:47
An identifier of a module&#39;s output.
Definition: InputTag_fwd.h:37
A class encapsulating a lua table.
Definition: ParameterSet.h:82
Definition: Pool.h:40
#define SQ(x)
Compute .
Definition: Math.h:25
A module declaration, defined from the configuration file.
Definition: Configuration.h:43