HepMC3 event record library
Feature.h
Go to the documentation of this file.
1// -*- C++ -*-
2//
3// This file is part of HepMC
4// Copyright (C) 2014-2019 The HepMC collaboration (see AUTHORS for details)
5//
6///
7/// @file Feature.h
8/// @brief Defines Feature interface for selecting Particles according to extracted Features.
9///
10
11#ifndef HEPMC3_FEATURE_H
12#define HEPMC3_FEATURE_H
13
14#include "HepMC3/GenParticle.h"
15#include "HepMC3/Filter.h"
16#include <functional>
17
18namespace HepMC3 {
19
20//////////////////////////////////////////////////////////////////////
21
22/**
23 * @brief GenericFeature defines the Feature interface
24 * GenericFeature is not intended to be used directly. The
25 * derived Feature class and its specialisations should be used.
26 *
27 * A Feature wraps a function object that can extract a
28 * generic Feature_type from a ConstGenParticlePtr. Usually the
29 * Feature_type would be something like int (e.g. status) or
30 * double (e.g. pT), but it could in principle be any attribute of a
31 * particle so long as there are well defined <, <=, >, >=, == and
32 * != operators for that attribute, as well as an abs function.
33 *
34 * Once a Feature is defined, you can obtain Filters that select
35 * Particles according to that Feature by e.g.
36 * Feature<int> status([](ConstGenParticlePtr p)->int{return p->status();});
37 * bool is_stable = (status == 1)(p);
38 * Filter is_beam = (status == 4);
39 * bool beam = is_beam(p);
40 *
41 * An abs function is also defined, so abs(Feature) works as you'd
42 * expect, e.g.
43 * Feature<double> rapidity([](ConstGenParticlePtr p)->double{return p->momentum().rap();});
44 * Filter rapCut = abs(rapidity) < 2.5;
45 *
46 * Please also see the Selector interface, which defines an
47 * abstract interface to Feature that is free of the template params
48 * and also includes some standard Features such as
49 *
50 * Selector::STATUS;
51 * Selector::PDG_ID;
52 * Selector::PT;
53 * Selector::RAPIDITY;
54 */
55template<typename Feature_type>
57
58public:
59
60 using Evaluator_type = std::function<Feature_type(ConstGenParticlePtr)>;
61 using EvaluatorPtr = std::shared_ptr<Evaluator_type>;
62
63 /// @brief access the underlying feature value
64 Feature_type operator()(ConstGenParticlePtr input)const {
65 return (*m_internal)(input);
66 }
67
68 /// @brief greater than operator
69 /// @return Filter function
70 Filter operator > (Feature_type value) const {
71 EvaluatorPtr functor = m_internal;
72 return [value, functor](ConstGenParticlePtr input)->bool{return (*functor)(input) > value;};
73 }
74 /// @brief less than operator
75 /// @return Filter function
76 Filter operator < (Feature_type value) const {
77 EvaluatorPtr functor = m_internal;
78 return [value, functor](ConstGenParticlePtr input)->bool{return (*functor)(input) < value;};
79 }
80
81 /// @brief greater than or equals operator
82 /// @return Filter function
83 Filter operator >= (Feature_type value) const {
84 EvaluatorPtr functor = m_internal;
85 return [value, functor](ConstGenParticlePtr input)->bool{return (*functor)(input) >= value;};
86 }
87
88 /// @brief less than or equals operator
89 /// @return Filter function
90 Filter operator <= (Feature_type value) const {
91 EvaluatorPtr functor = m_internal;
92 return [value, functor](ConstGenParticlePtr input)->bool{return (*functor)(input) <= value;};
93 }
94
95 /// @brief equality operator
96 /// @return Filter function
97 virtual Filter operator == (Feature_type value) const {
98 EvaluatorPtr functor = m_internal;
99 return [value, functor](ConstGenParticlePtr input)->bool{return (*functor)(input) == value;};
100 }
101
102 /// @brief inequality operator
103 /// @return Filter function
104 virtual Filter operator != (Feature_type value) const {
105 EvaluatorPtr functor = m_internal;
106 return [value, functor](ConstGenParticlePtr input)->bool{return (*functor)(input) != value;};
107 }
108
109protected:
110
111 /// Hide the constructor so no one can use GenericFeature directly
112 GenericFeature(Evaluator_type functor):m_internal(std::make_shared<Evaluator_type>(functor)) {}
113
114 /// Hide the copy constructor
115 GenericFeature(const GenericFeature &copy) : m_internal(copy.m_internal) {}
116
117 // internal copy of func for evaluation
118 // on the heap so will persist in resulting Filters even if
119 // parent Feature object was destroyed
120 EvaluatorPtr m_internal;
121
122};
123
124//////////////////////////////////////////////////////////////////////
125
126/** @brief Expose GenericFeature interface to derived Feature class
127 *
128 * This will get used for generic class types that aren't integral or
129 * floating point types.
130 *
131 * A Feature wraps a function object that can extract a
132 * generic Feature_type from a ConstGenParticlePtr. Usually the
133 * Feature_type would be something like int (e.g. status) or
134 * double (e.g. pT), but it could in principle be any attribute of a
135 * particle so long as there are well defined <, <=, >, >=, == and
136 * != operators for that attribute, as well as an abs function.
137 *
138 * Once a Feature is defined, you can obtain Filters that select
139 * Particles according to that Feature by e.g.
140 * Feature<int> status([](ConstGenParticlePtr p)->int{return p->status();});
141 * bool is_stable = (status == 1)(p);
142 * Filter is_beam = (status == 4);
143 * bool beam = is_beam(p);
144 *
145 * An abs function is also defined, so abs(Feature) works as you'd
146 * expect, e.g.
147 * Feature<double> rapidity([](ConstGenParticlePtr p)->double{return p->momentum().rap();});
148 * Filter rapCut = abs(rapidity) < 2.5;
149 *
150 * Please also see the Selector interface, which defines an
151 * abstract interface to Feature that is free of the template params
152 * and also includes some standard Features such as
153 *
154 * Selector::STATUS;
155 * Selector::PDG_ID;
156 * Selector::PT;
157 * Selector::RAPIDITY;
158
159 */
160template<typename Feature_type, typename Dummy=void>
161class Feature : public GenericFeature<Feature_type> {
162
163public:
164
165 using typename GenericFeature<Feature_type>::Evaluator_type;
166 using typename GenericFeature<Feature_type>::EvaluatorPtr;
167 using GenericFeature<Feature_type>::m_internal;
168
169 using GenericFeature<Feature_type>::operator ();
170 using GenericFeature<Feature_type>::operator >;
171 using GenericFeature<Feature_type>::operator >=;
172 using GenericFeature<Feature_type>::operator <;
173 using GenericFeature<Feature_type>::operator <=;
174 using GenericFeature<Feature_type>::operator ==;
175 using GenericFeature<Feature_type>::operator !=;
176
177 Feature(Evaluator_type functor) : GenericFeature<Feature_type>(functor) {}
178 Feature(const Feature &copy) : GenericFeature<Feature_type>(copy) {}
179
180 Feature<Feature_type> abs() const {
181 EvaluatorPtr functor = m_internal;
182 Evaluator_type absfunctor = [functor](ConstGenParticlePtr p)->Feature_type{return ::abs((*functor)(p));};
183 return Feature<Feature_type>(absfunctor);
184 }
185
186};
187
188//////////////////////////////////////////////////////////////////////
189
190/** @brief Specialisation of Feature for integral types
191 *
192 * It is a valid operator to compare an int to a float, but the
193 * generic version of these operators in the base class will
194 * first cast input float to an int, then compare that. In some cases
195 * the comparison will be incorrect because of rounding the float.
196 * e.g. int x=5; float y=5.5; bool result = x<y; would be wrong
197 * because y first gets converted to int 5.
198 *
199 * To solve this, we provide specialised comparison operators for
200 * integral type and double. Note that the opposite specialisation
201 * in which the Feature_type is floating_point is not necessary
202 */
203template<typename Feature_type>
204class Feature<Feature_type, typename std::enable_if<std::is_integral<Feature_type>::value, void>::type> : public GenericFeature<Feature_type> {
205
206public:
207
208 using GenericFeature<Feature_type>::operator ();
209 using GenericFeature<Feature_type>::operator >;
210 using GenericFeature<Feature_type>::operator >=;
211 using GenericFeature<Feature_type>::operator <;
212 using GenericFeature<Feature_type>::operator <=;
213 using GenericFeature<Feature_type>::operator ==;
214 using GenericFeature<Feature_type>::operator !=;
215
216 using typename GenericFeature<Feature_type>::Evaluator_type;
217 using typename GenericFeature<Feature_type>::EvaluatorPtr;
218
219 using GenericFeature<Feature_type>::m_internal;
220
221 Feature(Evaluator_type functor) : GenericFeature<Feature_type>(functor) {}
222 Feature(const Feature &copy) : GenericFeature<Feature_type>(copy) {}
223
224 Feature<Feature_type> abs() const {
225 EvaluatorPtr functor = m_internal;
226 Evaluator_type absfunctor = [functor](ConstGenParticlePtr p)->Feature_type{return ::abs((*functor)(p));};
227 return Feature<Feature_type>(absfunctor);
228 }
229
230 Filter operator > (double value) const {
231 EvaluatorPtr functor = m_internal;
232 return [value, functor](ConstGenParticlePtr input)->bool{return (*functor)(input) > value;};
233 }
234
235 Filter operator < (double value) const {
236 EvaluatorPtr functor = m_internal;
237 return [value, functor](ConstGenParticlePtr input)->bool{return (*functor)(input) < value;};
238 }
239
240 Filter operator == (double value) const {
241 EvaluatorPtr functor = m_internal;
242 return [value, functor](ConstGenParticlePtr input)->bool{
243 Feature_type local = (*functor)(input);
244 return fabs(local - value) <= ((::abs(local) < fabs(value))? fabs(value) : ::abs(local)) * std::numeric_limits<double>::epsilon();
245 };
246 }
247
248 Filter operator >= (double value) const { return !( (*this) < value );}
249
250 Filter operator <= (double value) const { return !( (*this) > value );}
251
252 Filter operator != (double value) const {
253 return !( (*this)==value );
254 }
255
256};
257
258//////////////////////////////////////////////////////////////////////
259
260/** @brief specialisation of Feature for floating point type
261 *
262 * Test of equality of floating point types is not safe. Here we
263 * provide a "reasonable" definition of equality based on the
264 * floating point precision.
265 */
266
267template<typename Feature_type>
268class Feature<Feature_type, typename std::enable_if<std::is_floating_point<Feature_type>::value, void>::type> : public GenericFeature<Feature_type> {
269
270public:
271
272 using typename GenericFeature<Feature_type>::Evaluator_type;
273 using typename GenericFeature<Feature_type>::EvaluatorPtr;
274
275 using GenericFeature<Feature_type>::operator ();
276 using GenericFeature<Feature_type>::operator >;
277 using GenericFeature<Feature_type>::operator >=;
278 using GenericFeature<Feature_type>::operator <;
279 using GenericFeature<Feature_type>::operator <=;
280
281 using GenericFeature<Feature_type>::m_internal;
282
283 Feature(Evaluator_type functor) : GenericFeature<Feature_type>(functor) {}
284 Feature(const Feature &copy) : GenericFeature<Feature_type>(copy) {}
285
286 Feature<Feature_type> abs() const {
287 EvaluatorPtr functor = m_internal;
288 Evaluator_type absfunctor = [functor](ConstGenParticlePtr p)->Feature_type{return fabs((*functor)(p));};
289 return Feature<Feature_type>(absfunctor);
290 }
291
292 Filter operator == (Feature_type value) const override {
293 EvaluatorPtr functor = m_internal;
294 return [value, functor](ConstGenParticlePtr input)->bool{
295 Feature_type local = (*functor)(input);
296 return fabs(local - value) <= ((fabs(local) < fabs(value))? fabs(value) : fabs(local)) * std::numeric_limits<Feature_type>::epsilon();
297 };
298 }
299
300 Filter operator != (Feature_type value) const override {
301 return !( (*this)==value );
302 }
303
304};
305
306//////////////////////////////////////////////////////////////////////
307
308/**
309 * @brief Obtain the absolute value of a Feature.
310 * This works as you'd expect. If foo is a valid Feature, then
311 * abs(foo) returns a new Feature that corresponds to the absolute
312 * value of the foo feature. You can construct a Filter from that in
313 * the usual way with e.g. Filter f = abs(foo) > 10.;
314 */
315template<typename Feature_type>
317 return input.abs();
318}
319
320}
321
322#endif
Defines Filter operations for combingin Filters.
Definition of class GenParticle.
Expose GenericFeature interface to derived Feature class.
Definition Feature.h:161
GenericFeature defines the Feature interface GenericFeature is not intended to be used directly....
Definition Feature.h:56
Filter operator>(Feature_type value) const
greater than operator
Definition Feature.h:70
Filter operator<(Feature_type value) const
less than operator
Definition Feature.h:76
Filter operator>=(Feature_type value) const
greater than or equals operator
Definition Feature.h:83
virtual Filter operator!=(Feature_type value) const
inequality operator
Definition Feature.h:104
GenericFeature(const GenericFeature &copy)
Hide the copy constructor.
Definition Feature.h:115
Feature_type operator()(ConstGenParticlePtr input) const
access the underlying feature value
Definition Feature.h:64
virtual Filter operator==(Feature_type value) const
equality operator
Definition Feature.h:97
Filter operator<=(Feature_type value) const
less than or equals operator
Definition Feature.h:90
GenericFeature(Evaluator_type functor)
Hide the constructor so no one can use GenericFeature directly.
Definition Feature.h:112
HepMC3 main namespace.
Definition ReaderGZ.h:28
std::function< bool(ConstGenParticlePtr)> Filter
type of Filter
Definition Filter.h:17
Feature< Feature_type > abs(const Feature< Feature_type > &input)
Obtain the absolute value of a Feature. This works as you'd expect. If foo is a valid Feature,...
Definition Feature.h:316