Belle II Software  release-05-01-25
threeHitFilters.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2014 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Jakob Lettenbichler (jakob.lettenbichler@oeaw.ac.at) *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <gtest/gtest.h>
12 
13 #include <tracking/spacePointCreation/SpacePoint.h>
14 #include <tracking/trackFindingVXD/filterMap/threeHitVariables/Angle3DSimple.h>
15 #include <tracking/trackFindingVXD/filterMap/threeHitVariables/Angle3DFull.h>
16 #include <tracking/trackFindingVXD/filterMap/threeHitVariables/CosAngleXY.h>
17 #include <tracking/trackFindingVXD/filterMap/threeHitVariables/AngleXYFull.h>
18 #include <tracking/trackFindingVXD/filterMap/threeHitVariables/AngleRZSimple.h>
19 #include <tracking/trackFindingVXD/filterMap/threeHitVariables/AngleRZFull.h>
20 #include <tracking/trackFindingVXD/filterMap/threeHitVariables/CircleDist2IP.h>
21 #include <tracking/trackFindingVXD/filterMap/threeHitVariables/SignCurvatureXY.h>
22 #include <tracking/trackFindingVXD/filterMap/threeHitVariables/SignCurvatureXYError.h>
23 #include <tracking/trackFindingVXD/filterMap/threeHitVariables/DeltaSlopeRZ.h>
24 #include <tracking/trackFindingVXD/filterMap/threeHitVariables/DeltaSlopeZoverS.h>
25 #include <tracking/trackFindingVXD/filterMap/threeHitVariables/DeltaSoverZ.h>
26 #include <tracking/trackFindingVXD/filterMap/threeHitVariables/HelixParameterFit.h>
27 #include <tracking/trackFindingVXD/filterMap/threeHitVariables/Pt.h>
28 #include <tracking/trackFindingVXD/filterMap/threeHitVariables/CircleRadius.h>
29 #include <tracking/trackFindingVXD/filterTools/SelectionVariableHelper.h>
30 
31 #include <tracking/trackFindingVXD/filterMap/filterFramework/Shortcuts.h>
32 
33 #include <vxd/geometry/SensorInfoBase.h>
34 
35 #include <math.h>
36 
37 
38 using namespace std;
39 using namespace Belle2;
40 
41 namespace VXDTFthreeHitFilterTest {
42 
44  class ThreeHitFilterTest : public ::testing::Test {
45  protected:
46  };
47 
48 
49 
51  VXD::SensorInfoBase createSensorInfo(VxdID aVxdID, double globalX = 0., double globalY = 0., double globalZ = -0.)
52  {
53  // (SensorType type, VxdID id, double width, double length, double thickness, int uCells, int vCells, double width2=-1, double splitLength=-1, int vCells2=0)
54  VXD::SensorInfoBase sensorInfoBase(VXD::SensorInfoBase::PXD, aVxdID, 2.3, 4.2, 0.3, 2, 4, -1);
55 
56  TGeoRotation r1;
57  r1.SetAngles(45, 20, 30); // rotation defined by Euler angles
58  TGeoTranslation t1(globalX, globalY, globalZ);
59  TGeoCombiTrans c1(t1, r1);
60  TGeoHMatrix transform = c1;
61  sensorInfoBase.setTransformation(transform);
62  // also need reco-transform
63  sensorInfoBase.setTransformation(transform, true);
64 
65  return sensorInfoBase;
66  }
67 
68 
69 
71  SpacePoint provideSpacePointDummy(double X, double Y, double Z)
72  {
73  VxdID aVxdID = VxdID(1, 1, 1);
74  VXD::SensorInfoBase sensorInfoBase = createSensorInfo(aVxdID, X, Y, Z);
75 
76  PXDCluster aCluster = PXDCluster(aVxdID, 0., 0., 0.1, 0.1, 0, 0, 1, 1, 1, 1, 1, 1);
77 
78  return SpacePoint(&aCluster, &sensorInfoBase);
79  }
80 
81 
83  SpacePoint provideSpacePointDummyError(double X, double Y, double Z, double eX, double eY, double eZ)
84  {
85  VxdID aVxdID = VxdID(1, 1, 1);
86 
87  B2Vector3D pos(X, Y, Z);
88  B2Vector3D poserr(eX, eY, eZ);
89 
90  return SpacePoint(pos, poserr, {0.0, 0.0}, {false, false}, aVxdID, Belle2::VXD::SensorInfoBase::PXD);
91  }
92 
93 
94 
96  double lastResult = 0.;
97 
98 
99 
101  class ResultsObserver : public VoidObserver {
102  public:
104  template<class Var, typename ... otherTypes>
105  static void notify(const Var& filterType,
106  typename Var::variableType fResult,
107  otherTypes ...)
108  {
109  B2INFO("ResultsObserver: Filter " << filterType.name() << " got result of " << fResult);
110  lastResult = fResult;
111  }
112 
113  };
114 
115 
117  TEST_F(ThreeHitFilterTest, SpacePointCreation)
118  {
119  SpacePoint testSP = provideSpacePointDummy(1.2, 2.3, 4.2);
120  EXPECT_FLOAT_EQ(1.2, testSP.getPosition().X()) ;
121  EXPECT_FLOAT_EQ(2.3, testSP.getPosition().Y()) ;
122  EXPECT_FLOAT_EQ(4.2, testSP.getPosition().Z()) ;
123  }
124 
125 
127  TEST_F(ThreeHitFilterTest, SelectionVariableName)
128  {
129  auto angle3DS = Angle3DSimple<SpacePoint>();
130  EXPECT_EQ("Angle3DSimple" , angle3DS.name());
131  auto angle3DF = Angle3DFull<SpacePoint>();
132  EXPECT_EQ("Angle3DFull" , angle3DF.name());
133  auto angleXYF = AngleXYFull<SpacePoint>();
134  EXPECT_EQ("AngleXYFull" , angleXYF.name());
135  auto angleXYS = CosAngleXY<SpacePoint>();
136  EXPECT_EQ("CosAngleXY" , angleXYS.name());
137  auto angleRZF = AngleRZFull<SpacePoint>();
138  EXPECT_EQ("AngleRZFull" , angleRZF.name());
139  auto angleRZS = AngleRZSimple<SpacePoint>();
140  EXPECT_EQ("AngleRZSimple" , angleRZS.name());
141  auto dSlopeRZ = DeltaSlopeRZ<SpacePoint>();
142  EXPECT_EQ("DeltaSlopeRZ" , dSlopeRZ.name());
143  auto helixFit = HelixParameterFit<SpacePoint>();
144  EXPECT_EQ("HelixParameterFit" , helixFit.name());
145  auto dSlopeZS = DeltaSlopeZoverS<SpacePoint>();
146  EXPECT_EQ("DeltaSlopeZoverS" , dSlopeZS.name());
147  auto dSZ = DeltaSoverZ<SpacePoint>();
148  EXPECT_EQ("DeltaSoverZ" , dSZ.name());
149  auto pT = Pt<SpacePoint>();
150  EXPECT_EQ("Pt" , pT.name());
151  auto cDist2IP = CircleDist2IP<SpacePoint>();
152  EXPECT_EQ("CircleDist2IP" , cDist2IP.name());
153  }
154 
155 
157  TEST_F(ThreeHitFilterTest, BasicFilterTestAngle3DSimple)
158  {
159  // Very verbose declaration, see below for convenient shortcuts
161 
162  SpacePoint x1 = provideSpacePointDummy(0. , 0., 0.);
163  SpacePoint x2 = provideSpacePointDummy(.5 , 0., 0.);
164  SpacePoint x3 = provideSpacePointDummy(2. , 0., 0.);
165 
166  EXPECT_FALSE(filter.accept(x1, x2, x3));
167 
168  }
169 
170 
172  TEST_F(ThreeHitFilterTest, TestAngles)
173  {
174  SpacePoint outerSP = provideSpacePointDummy(6., 4., 1.);
175  SpacePoint centerSP = provideSpacePointDummy(3., 3., 0.);
176  SpacePoint innerSP = provideSpacePointDummy(1., 1., 0.);
177 
178 
180  0.091));
181  EXPECT_TRUE(filterAngle3D.accept(outerSP, centerSP, innerSP));
182  EXPECT_FLOAT_EQ(0.090909090909090909091, lastResult); // last result is what filterAngle3D had calculated
183 
184  Filter< Angle3DFull<SpacePoint>, Range<double, double>, ResultsObserver > filterAngle3Dfull(Range<double, double>(31.48,
185  31.49));
186  EXPECT_TRUE(filterAngle3Dfull.accept(outerSP, centerSP, innerSP));
187  EXPECT_FLOAT_EQ(31.4821541052938556040832384555411729852856, lastResult);
188 
189 
190  Filter< AngleXYFull<SpacePoint>, Range<double, double>, ResultsObserver > filterAngleXYfull(Range<double, double>(26.5,
191  26.6));
192  EXPECT_TRUE(filterAngleXYfull.accept(outerSP, centerSP, innerSP));
193  EXPECT_FLOAT_EQ(26.5650511770779893515721937204532946712042, lastResult);
194 
195  Filter< CosAngleXY<SpacePoint>, Range<double, double>, ResultsObserver > filterAngleXY(Range<double, double>(0.7,
196  0.9));
197  EXPECT_TRUE(filterAngleXY.accept(outerSP, centerSP, innerSP));
198  EXPECT_NEAR(0.89442718, lastResult, 0.0000001);
199 
200 
201  Filter< AngleRZFull<SpacePoint>, Range<double, double>, ResultsObserver > filterAngleRZfull(Range<double, double>(17.5,
202  17.6));
203  EXPECT_TRUE(filterAngleRZfull.accept(outerSP, centerSP, innerSP));
204  EXPECT_NEAR(17.548400613792298064, lastResult, 0.001);
205 
207  0.96));
208  EXPECT_TRUE(filterAngleRZ.accept(outerSP, centerSP, innerSP));
209  EXPECT_NEAR(cos(17.54840061379229806435203716652846677620 * M_PI / 180.), lastResult, 0.001);
210  }
211 
212 
213 
215  TEST_F(ThreeHitFilterTest, TestSignAndOtherFilters)
216  {
217  SpacePoint outerSP = provideSpacePointDummy(6., 4., 1.);
218  SpacePoint centerSP = provideSpacePointDummy(3., 3., 0.);
219  SpacePoint innerSP = provideSpacePointDummy(1., 1., 0.);
220  SpacePoint outerHighSP = provideSpacePointDummy(4., 6., 1.);
221  B2Vector3D sigma(.01, .01, .01), unrealsigma(2, 2, 2);
222 
223  SpacePoint outerSPsigma = provideSpacePointDummyError(6., 4., 1., 0.01, 0.01, 0.01);
224  SpacePoint centerSPsigma = provideSpacePointDummyError(3., 3., 0., 0.01, 0.01, 0.01);
225  SpacePoint innerSPsigma = provideSpacePointDummyError(1., 1., 0., 0.01, 0.01, 0.01);
226  SpacePoint outerHighSPsigma = provideSpacePointDummyError(4., 6., 1., 0.01, 0.01, 0.01);
227 
228  SpacePoint outerSPunrealsigma = provideSpacePointDummyError(6., 4., 1., 2., 2., 2.);
229  SpacePoint centerSPunrealsigma = provideSpacePointDummyError(3., 3., 0., 2., 2., 2.);
230  SpacePoint innerSPunrealsigma = provideSpacePointDummyError(1., 1., 0., 2., 2., 2.);
231 
232  Filter< DeltaSlopeRZ<SpacePoint>, Range<double, double>, ResultsObserver > filterDeltaSlopeRZ(Range<double, double>(0.3,
233  0.31));
234  EXPECT_TRUE(filterDeltaSlopeRZ.accept(outerSP, centerSP, innerSP));
235  EXPECT_FLOAT_EQ(0.30627736916966945608, lastResult);
236 
238  0.01));
239  EXPECT_TRUE(filterHelixFit.accept(outerSP, centerSP, innerSP));
240  EXPECT_FLOAT_EQ(0., lastResult);
241 
242  EXPECT_DOUBLE_EQ(1., SignCurvatureXYError<SpacePoint>::value(outerSPsigma, centerSPsigma, innerSPsigma));
243  EXPECT_DOUBLE_EQ(-1., SignCurvatureXYError<SpacePoint>::value(outerHighSPsigma, centerSPsigma, innerSPsigma));
244  EXPECT_DOUBLE_EQ(-1., SignCurvatureXYError<SpacePoint>::value(innerSPsigma, centerSPsigma, outerSPsigma));
245  //for very large sigma, this track is approximately straight:
246  EXPECT_DOUBLE_EQ(0., SignCurvatureXYError<SpacePoint>::value(outerSPunrealsigma, centerSPunrealsigma, innerSPunrealsigma));
247 
248  EXPECT_LT(0., SignCurvatureXY<SpacePoint>::value(outerSP, centerSP, innerSP));
249  EXPECT_GT(0., SignCurvatureXY<SpacePoint>::value(outerHighSP, centerSP, innerSP));
250  EXPECT_GT(0., SignCurvatureXY<SpacePoint>::value(innerSP, centerSP, outerSP));
251  EXPECT_LT(0., SignCurvatureXY<SpacePoint>::value(outerSP, centerSP, innerSP));
252 
253  EXPECT_DOUBLE_EQ(1., SignCurvatureXY<SpacePoint>::value(outerSP, centerSP, innerSP));
254  EXPECT_DOUBLE_EQ(-1., SignCurvatureXY<SpacePoint>::value(outerHighSP, centerSP, innerSP));
255  EXPECT_DOUBLE_EQ(-1., SignCurvatureXY<SpacePoint>::value(innerSP, centerSP, outerSP));
256  }
257 
258 
259 
261  TEST_F(ThreeHitFilterTest, TestDeltaSOverZ)
262  {
263  lastResult = 0;
264  B2Vector3D iVec(1., 1., 0.), cVec(3., 3., 1.), oVec(6., 4., 3.);
265  SpacePoint outerSP = provideSpacePointDummy(oVec.X(), oVec.Y(), oVec.Z());
266  SpacePoint centerSP = provideSpacePointDummy(cVec.X(), cVec.Y(), cVec.Z());
267  SpacePoint innerSP = provideSpacePointDummy(iVec.X(), iVec.Y(), iVec.Z());
268 
269  Filter< DeltaSlopeZoverS<SpacePoint>, Range<double, double>, ResultsObserver > filterDeltaSlopeZOverS(Range<double, double>
270  (0.31,
271  0.32));
272  EXPECT_TRUE(filterDeltaSlopeZOverS.accept(outerSP, centerSP, innerSP));
273  EXPECT_NEAR(0.31823963, lastResult, 0.00001);
274 
275  Filter< DeltaSoverZ<SpacePoint>, Range<double, double>, ResultsObserver > filterDeltaSOverZ(Range<double, double>(-0.39,
276  -0.38));
277  EXPECT_TRUE(filterDeltaSOverZ.accept(outerSP, centerSP, innerSP));
278  EXPECT_FLOAT_EQ(-0.38471845, lastResult);
279 
280 
281  //calcDeltaSOverZV2 is invariant under rotations in the r-z plane:
282  oVec.RotateZ(.4);
283  cVec.RotateZ(.4);
284  iVec.RotateZ(.4);
285  SpacePoint outerSP2 = provideSpacePointDummy(oVec.X(), oVec.Y(), oVec.Z());
286  SpacePoint centerSP2 = provideSpacePointDummy(cVec.X(), cVec.Y(), cVec.Z());
287  SpacePoint innerSP2 = provideSpacePointDummy(iVec.X(), iVec.Y(), iVec.Z());
288 
289  EXPECT_TRUE(filterDeltaSlopeZOverS.accept(outerSP2, centerSP2, innerSP2));
290  EXPECT_NEAR(0.31823963, lastResult, 0.00001);
291 
292  EXPECT_TRUE(filterDeltaSOverZ.accept(outerSP2, centerSP2, innerSP2));
293  EXPECT_FLOAT_EQ(-0.38471845, lastResult);
294  }
295 
296 
297 
299  TEST_F(ThreeHitFilterTest, TestCalcPt)
300  {
301  typedef SelVarHelper<SpacePoint, double> Helper;
302 
303  SpacePoint outerSP = provideSpacePointDummy(3., 4., 3.);
304  SpacePoint centerSP = provideSpacePointDummy(3., 2., 1.);
305  SpacePoint innerSP = provideSpacePointDummy(1., 2., 0.);
306 
307  SpacePoint outerSPEvil = provideSpacePointDummy(6., 3., 0.);
308  SpacePoint centerSP2 = provideSpacePointDummy(3., 3., 0.);
309  SpacePoint innerSP2 = provideSpacePointDummy(1., 1., 0.);
310 
311  double ptTrue = 0,
312  rawRadius = 1.414213562373095048801688724209698078570;
313  ptTrue = Helper::calcPt(rawRadius);
314  Filter< Pt<SpacePoint>, Range<double, double>, ResultsObserver > filterPt(Range<double, double>(0.005, 0.05));
315  EXPECT_TRUE(filterPt.accept(outerSP, centerSP, innerSP));
316  EXPECT_FLOAT_EQ(ptTrue, lastResult);
317 
318  ptTrue = 0.017118925181688543;
319  EXPECT_TRUE(filterPt.accept(outerSPEvil, centerSP2, innerSP2));
320  EXPECT_FLOAT_EQ(ptTrue, lastResult);
321 
322  //B2WARNING("MUST produce errors: 2 hits are the same: " << ptTrue << ", Pt: " << pt );
323  EXPECT_ANY_THROW(filterPt.accept(outerSP, outerSP, innerSP));
324 
325  Helper::resetMagneticField(0.976);
326  double ptTrue1p5 = ptTrue;
327  ptTrue = Helper::calcPt(rawRadius);
328  EXPECT_FALSE(filterPt.accept(outerSP, centerSP, innerSP));
329  EXPECT_FLOAT_EQ(ptTrue, lastResult);
330  EXPECT_LT(ptTrue, ptTrue1p5);
331 
332  Helper::resetMagneticField(1.5);
333  }
334 
335 
337  TEST_F(ThreeHitFilterTest, TestCircleRadius)
338  {
339  SpacePoint outerSP = provideSpacePointDummy(3., 4., 3.);
340  SpacePoint centerSP = provideSpacePointDummy(3., 2., 1.);
341  SpacePoint innerSP = provideSpacePointDummy(1., 2., 0.);
342 
343  double rawRadius = 1.414213562373095048801688724209698078570;
344  Filter< CircleRadius<SpacePoint>, Range<double, double>, ResultsObserver > filterCrad(Range<double, double>(1.4, 1.42));
345  EXPECT_TRUE(filterCrad.accept(outerSP, centerSP, innerSP));
346  EXPECT_FLOAT_EQ(rawRadius, lastResult);
347 
348  //B2WARNING("MUST produce errors: 2 hits are the same:");
349  EXPECT_ANY_THROW(filterCrad.accept(outerSP, outerSP, innerSP));
350 
351 // EXPECT_FALSE(filterCrad.accept(outerSP, centerSP, innerSP));
352  EXPECT_FLOAT_EQ(rawRadius, lastResult); // after fail result is not changed
353  }
354 
356  TEST_F(ThreeHitFilterTest, TestCircleDist2IP)
357  {
358  SpacePoint outerSPEvil = provideSpacePointDummy(6., 3., 0.);
359  SpacePoint outerSPSimple = provideSpacePointDummy(6., 4., 0.);
360  SpacePoint centerSP = provideSpacePointDummy(3., 3., 0.);
361  SpacePoint innerSP = provideSpacePointDummy(1., 1., 0.);
362 
363  Filter< CircleDist2IP<SpacePoint>, Range<double, double>, ResultsObserver > filteCircleDist2IP(Range<double, double>(0.44,
364  0.45));
365  EXPECT_TRUE(filteCircleDist2IP.accept(outerSPSimple, centerSP, innerSP));
366  EXPECT_NEAR(0.44499173338941, lastResult, 0.00001);
367 
368  EXPECT_FALSE(filteCircleDist2IP.accept(outerSPEvil, centerSP, innerSP));
369  EXPECT_FLOAT_EQ(0.719806016136754, lastResult);
370 
371 
372  typedef SelVarHelper<SpacePoint, double> Helper;
373  Helper::resetMagneticField(0.976);
374 
375  EXPECT_TRUE(filteCircleDist2IP.accept(outerSPSimple, centerSP, innerSP));
376  EXPECT_NEAR(0.44499173338941, lastResult, 0.00001);
377 
378  EXPECT_FALSE(filteCircleDist2IP.accept(outerSPEvil, centerSP, innerSP));
379  EXPECT_FLOAT_EQ(0.719806016136754, lastResult);
380 
381  Helper::resetMagneticField(1.5);
382  }
383 
384 }
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
Belle2::SelVarHelper
contains a collection of functions and related stuff needed for SelectionVariables implementing 2-,...
Definition: SelectionVariableHelper.h:31
Belle2::VXD::SensorInfoBase
Base class to provide Sensor Information for PXD and SVD.
Definition: SensorInfoBase.h:40
Belle2::filter
std::map< ExpRun, std::pair< double, double > > filter(const std::map< ExpRun, std::pair< double, double >> &runs, double cut, std::map< ExpRun, std::pair< double, double >> &runsRemoved)
filter events to remove runs shorter than cut, it stores removed runs in runsRemoved
Definition: Splitter.cc:43
Belle2::B2Vector3::Z
DataType Z() const
access variable Z (= .at(2) without boundary check)
Definition: B2Vector3.h:434
Belle2::SpacePoint
SpacePoint typically is build from 1 PXDCluster or 1-2 SVDClusters.
Definition: SpacePoint.h:52
Belle2::SpacePoint::getPosition
const B2Vector3< double > & getPosition() const
return the position vector in global coordinates
Definition: SpacePoint.h:148
Belle2::B2Vector3< double >
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::Range
Represents a range of arithmetic types.
Definition: Range.h:39
Belle2::TEST_F
TEST_F(GlobalLabelTest, LargeNumberOfTimeDependentParameters)
Test large number of time-dep params for registration and retrieval.
Definition: globalLabel.cc:65
VXDTFthreeHitFilterTest::ThreeHitFilterTest
Test class for these new and shiny two-hit-filters.
Definition: threeHitFilters.cc:44
Belle2::PXDCluster
The PXD Cluster class This class stores all information about reconstructed PXD clusters The position...
Definition: PXDCluster.h:41
VXDTFthreeHitFilterTest::ResultsObserver
takes result, prints it and stores it to lastResult
Definition: threeHitFilters.cc:101
Belle2::VXD::SensorInfoBase::PXD
@ PXD
PXD Sensor.
Definition: SensorInfoBase.h:44
Belle2::B2Vector3::X
DataType X() const
access variable X (= .at(0) without boundary check)
Definition: B2Vector3.h:430
VXDTFthreeHitFilterTest::ResultsObserver::notify
static void notify(const Var &filterType, typename Var::variableType fResult, otherTypes ...)
notify function is called by the filter, this one takes result, prints it and stores it to lastResult...
Definition: threeHitFilters.cc:105
Belle2::Filter
This class is used to select pairs, triplets...
Definition: Filter.h:44
Belle2::B2Vector3::Y
DataType Y() const
access variable Y (= .at(1) without boundary check)
Definition: B2Vector3.h:432
Belle2::VoidObserver
The most CPU efficient Observer for the VXDTF filter tools (even if useless).
Definition: VoidObserver.h:40