Jpp
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Public Types | Public Member Functions | Public Attributes | Static Public Attributes | List of all members
JFIT::JRegressor< JShower3EZ, JGandalf > Struct Template Reference

Regressor function object for JShower3EZ fit using JGandalf minimiser. More...

#include <JShower3EZRegressor.hh>

Inheritance diagram for JFIT::JRegressor< JShower3EZ, JGandalf >:
JFIT::JAbstractRegressor< JShower3EZ, JGandalf > JFIT::JGandalf< JShower3EZ > JEEP::JMessage< T >

Public Types

typedef JTOOLS::JSplineFunction1S_t JFunction1D_t
 
typedef JTOOLS::JMapList
< JTOOLS::JPolint0FunctionalMap,
JTOOLS::JMapList
< JTOOLS::JPolint0FunctionalMap,
JTOOLS::JMapList
< JTOOLS::JPolint0FunctionalGridMap,
JTOOLS::JMapList
< JTOOLS::JPolint0FunctionalGridMap > > > > 
JPDFMaplist_t
 
typedef JPHYSICS::JPDFTable
< JFunction1D_t, JPDFMaplist_t
JPDF_t
 
typedef JTOOLS::JMapList
< JTOOLS::JPolint1FunctionalMapH,
JTOOLS::JMapList
< JTOOLS::JPolint1FunctionalMap,
JTOOLS::JMapList
< JTOOLS::JPolint1FunctionalGridMap,
JTOOLS::JMapList
< JTOOLS::JPolint1FunctionalGridMap > > > > 
JNPEMaplist_t
 
typedef JPHYSICS::JNPETable
< double, double,
JNPEMaplist_t
JNPE_t
 
typedef JGandalf< JShower3EZminimiser_type
 
typedef JRegressor< JShower3EZ,
JGandalf
regressor_type
 
typedef minimiser_type::result_type result_type
 
typedef JShower3EZ::parameter_type parameter_type
 Data type of fit parameter. More...
 

Public Member Functions

 JRegressor (const std::string &fileDescriptor)
 Parameterized constructor. More...
 
result_type operator() (const JShower3EZ &shower, const JPMTW0 &pmt) const
 Fit function. More...
 
JNPE_t::result_type getH0 (const double R_Hz) const
 Get background hypothesis value for time integrated PDF. More...
 
JNPE_t::result_type getH1 (const double D, const double cosDelta, const double theta, const double phi, const double E) const
 Get signal hypothesis value for time integrated PDF. More...
 
double getRmax () const
 Get maximal road width of PDF. More...
 
result_type operator() (const JShower3EZ &value, T __begin, T __end)
 Global fit. More...
 
result_type operator() (const JFunction_t &fit, T __begin, T __end, Args...args)
 Multi-dimensional fit of two data sets. More...
 

Public Attributes

JNPE_t npe [NUMBER_OF_PDFS]
 PDF. More...
 
JLANG::JSharedPointer
< JMEstimator
estimator
 M-Estimator function. More...
 
double lambda
 
JShower3EZ value
 
JShower3EZ error
 
std::vector< parameter_typeparameters
 
int numberOfIterations
 
JMATH::JMatrixNS H
 

Static Public Attributes

static JTimeRange T_ns
 Time window with respect to Cherenkov hypothesis [ns]. More...
 
static double Vmax_npe = std::numeric_limits<double>::max()
 Maximal integral of PDF [npe]. More...
 
static const int NUMBER_OF_PDFS = 2
 
static const JFIT::JPDFType_t pdf_t [NUMBER_OF_PDFS]
 
static int MAXIMUM_ITERATIONS
 maximal number of iterations More...
 
static double EPSILON
 maximal distance to minimum More...
 
static double LAMBDA_MIN
 minimal value control parameter More...
 
static double LAMBDA_MAX
 maximal value control parameter More...
 
static double LAMBDA_UP
 multiplication factor control parameter More...
 
static double LAMBDA_DOWN
 multiplication factor control parameter More...
 
static double PIVOT
 minimal value diagonal element of matrix More...
 
static int debug = 0
 debug level (default is off). More...
 

Detailed Description

template<>
struct JFIT::JRegressor< JShower3EZ, JGandalf >

Regressor function object for JShower3EZ fit using JGandalf minimiser.

Definition at line 253 of file JShower3EZRegressor.hh.

Member Typedef Documentation

Definition at line 258 of file JShower3EZRegressor.hh.

Definition at line 262 of file JShower3EZRegressor.hh.

Definition at line 263 of file JShower3EZRegressor.hh.

Definition at line 268 of file JShower3EZRegressor.hh.

Definition at line 269 of file JShower3EZRegressor.hh.

Definition at line 78 of file JRegressor.hh.

Definition at line 79 of file JRegressor.hh.

typedef minimiser_type::result_type JFIT::JAbstractRegressor< JShower3EZ , JGandalf >::result_type
inherited

Definition at line 80 of file JRegressor.hh.

Data type of fit parameter.

Definition at line 63 of file JGandalf.hh.

Constructor & Destructor Documentation

JFIT::JRegressor< JShower3EZ, JGandalf >::JRegressor ( const std::string &  fileDescriptor)
inline

Parameterized constructor.

The PDF file descriptor should contain the wild card character JPHYSICS::WILD_CARD which will be replaced by the corresponding PDF types listed in JRegressor<JShower3Z, JGandalf>::pdf_t.

Parameters
fileDescriptorPDF file descriptor

Definition at line 280 of file JShower3EZRegressor.hh.

280  :
281  estimator(new JMEstimatorNull())
282  {
283  using namespace std;
284  using namespace JPP;
285 
286  const JPDF_t::JSupervisor supervisor(new JPDF_t::JDefaultResult(JMATH::zero));
287 
288  for (int i = 0; i != NUMBER_OF_PDFS; ++i) {
289 
290  try {
291 
292  JPDF_t pdf;
293 
294  const string file_name = getFilename(fileDescriptor, pdf_t[i]);
295 
296  NOTICE("loading PDF from file " << file_name << "... " << flush);
297 
298  pdf.load(file_name.c_str());
299 
300  NOTICE("OK" << endl);
301 
302  pdf.setExceptionHandler(supervisor);
303 
304  npe[i] = JNPE_t(pdf);
305  }
306  catch(const JException& error) {
307  FATAL(error.what() << endl);
308  }
309  }
310 
311  // Add PDFs
312  for (int i = 1; i < NUMBER_OF_PDFS; i += 2) {
313 
314  npe[ i ].add(npe[i-1]);
315 
316  JNPE_t buffer;
317 
318  npe[i-1].swap(buffer);
319  }
320  }
JPHYSICS::JPDFTable< JFunction1D_t, JPDFMaplist_t > JPDF_t
static const JZero zero
Function object to assign zero value.
Definition: JZero.hh:105
static const JFIT::JPDFType_t pdf_t[NUMBER_OF_PDFS]
void add(const JNPETable &input)
Add NPE table.
Definition: JNPETable.hh:124
JPHYSICS::JNPETable< double, double, JNPEMaplist_t > JNPE_t
#define NOTICE(A)
Definition: JMessage.hh:64
#define FATAL(A)
Definition: JMessage.hh:67
std::string getFilename(const std::string &file_name)
Get file name part, i.e. part after last JEEP::PATHNAME_SEPARATOR if any.
Definition: JeepToolkit.hh:88
JLANG::JSharedPointer< JMEstimator > estimator
M-Estimator function.

Member Function Documentation

result_type JFIT::JRegressor< JShower3EZ, JGandalf >::operator() ( const JShower3EZ shower,
const JPMTW0 pmt 
) const
inline

Fit function.

This method is used to determine the chi2 of given PMT with respect to shower hypothesis.

Parameters
showershower
pmtpmt
Returns
chi2

Definition at line 330 of file JShower3EZRegressor.hh.

331  {
332  using namespace JPP;
333  using namespace std;
334 
335  JPosition3D D(pmt.getPosition());
336  JDirection3D U(pmt.getDirection());
337 
338  D.sub(shower.getPosition());
339 
340  JVersor3D shower_dir(shower.getDirection());
341 
342  const double z = D.getDot(shower_dir);
343  const double x = D.getX() - z * shower.getDX();
344  const double y = D.getY() - z * shower.getDY();
345  const double R2 = D.getLengthSquared() - z*z;
346  const double R = sqrt(R2);
347  const double d = D.getLength();
348  const double cos_th0 = z/d; // theta0 = angle between shower direction and PMT position
349 
350  U.rotate(JRotation3Z(-atan2(y,x))); // rotate PMT axis to x-z plane
351 
352  const double theta = U.getTheta();
353  const double phi = fabs(U.getPhi());
354 
355  JNPE_t::result_type H0 = getH0(pmt.getR()); // getH0 = Get background hypothesis value for time integrated PDF.
356  JNPE_t::result_type H1 = getH1(d, cos_th0, theta, phi, shower.getE()); // getH1 = Get signal hypothesis value for time integrated PDF.
357 
358  if (get_value(H1) >= Vmax_npe) {
359  H1 *= Vmax_npe / get_value(H1);
360  }
361 
362  H1 += H0; // now H1 is signal + background
363 
364  const bool hit = pmt.getN() != 0;
365  const double u = H1.getChi2(hit);
366 
368 
369  result.chi2 = estimator->getRho(u);
370 
371  result.gradient = JShower3EZ(JPoint4D(JVector3D(x*z / (d*d*d), // d(cos_th0)/d(x)
372  y*z / (d*d*d), // d(cos_th0)/d(y)
373  -(R*R) / (d*d*d)), // d(cos_th0)/d(z)
374  0.0), // d(cos_th0)/d(t)
375  JVersor3Z(x*z*z / (d*R), // d(cos_th0)/d(dx)
376  y*z*z / (d*R)), // d(cos_th0)/d(dy)
377  0.0); // d(cos_th0)/d(E)
378 
379  result.gradient.mul(estimator->getPsi(u));
380  result.gradient.mul(H1.getDerivativeOfChi2(hit)); // x d(chi2)/d(cos_th0)
381 
382  return result;
383  }
do echo Generating $dir eval D
Definition: JDrawLED.sh:50
esac $JPP_DIR examples JDetector JTransitTime o $OUTPUT_FILE n N $NPE T $TTS_NS d $DEBUG for HISTOGRAM in tts tt2 pmt
Definition: JTransitTime.sh:36
return result
Definition: JPolint.hh:695
then print_variable DETECTOR INPUT_FILE INTERMEDIATE_FILE check_input_file $DETECTOR $INPUT_FILE check_output_file $INTERMEDIATE_FILE $OUTPUT_FILE JMCEvt f $INPUT_FILE o $INTERMEDIATE_FILE d
Definition: JPath.sh:52
then usage $script[distance] fi case set_variable R
Definition: JDrawLED.sh:40
JNPE_t::result_type getH1(const double D, const double cosDelta, const double theta, const double phi, const double E) const
Get signal hypothesis value for time integrated PDF.
static double Vmax_npe
Maximal integral of PDF [npe].
JLANG::JSharedPointer< JMEstimator > estimator
M-Estimator function.
double u[N+1]
Definition: JPolint.hh:706
JResultEvaluator< JResult_t >::result_type get_value(const JResult_t &value)
Helper method to recursively evaluate a to function value.
Definition: JResult.hh:936
JNPE_t::result_type getH0(const double R_Hz) const
Get background hypothesis value for time integrated PDF.
JNPE_t::result_type JFIT::JRegressor< JShower3EZ, JGandalf >::getH0 ( const double  R_Hz) const
inline

Get background hypothesis value for time integrated PDF.

Parameters
R_Hzrate [Hz]
Returns
hypothesis value

Definition at line 391 of file JShower3EZRegressor.hh.

392  {
393  return JNPE_t::result_type(R_Hz * 1e-9 * T_ns.getLength(), 0.0);
394  }
static JTimeRange T_ns
Time window with respect to Cherenkov hypothesis [ns].
JNPE_t::result_type JFIT::JRegressor< JShower3EZ, JGandalf >::getH1 ( const double  D,
const double  cosDelta,
const double  theta,
const double  phi,
const double  E 
) const
inline

Get signal hypothesis value for time integrated PDF.

Parameters
DPMT distance from shower [m]
cosDeltaangle between shower direction and PMT position
thetaPMT zenith angle [deg]
phiPMT azimuth angle [deg]
Eshower energy [GeV]
Returns
hypothesis value

Definition at line 406 of file JShower3EZRegressor.hh.

411  {
412  JNPE_t::result_type h1 = JMATH::zero;
413 
414  for (int i = 0; i != NUMBER_OF_PDFS; ++i) {
415 
416  if (!npe[i].empty() && D <= npe[i].getXmax()) {
417 
418  try {
419 
420  JNPE_t::result_type y1 = E * npe[i](std::max(D, npe[i].getXmin()), cosDelta, theta, phi);
421 
422  if (get_value(y1) > 0.0) {
423  h1 += y1;
424  }
425 
426  }
427  catch(JLANG::JException& error) {
428  ERROR(error << std::endl);
429  }
430  }
431  }
432 
433  return h1;
434  }
General exception.
Definition: JException.hh:23
do echo Generating $dir eval D
Definition: JDrawLED.sh:50
then for HISTOGRAM in h0 h1
Definition: JMatrixNZ.sh:69
static const JZero zero
Function object to assign zero value.
Definition: JZero.hh:105
JResultEvaluator< JResult_t >::result_type get_value(const JResult_t &value)
Helper method to recursively evaluate a to function value.
Definition: JResult.hh:936
then usage $script[input file[working directory[option]]] nWhere option can be E
Definition: JMuonPostfit.sh:37
double JFIT::JRegressor< JShower3EZ, JGandalf >::getRmax ( ) const
inline

Get maximal road width of PDF.

Returns
road width [m]

Definition at line 442 of file JShower3EZRegressor.hh.

443  {
444  double xmax = 0.0;
445 
446  for (int i = 0; i != NUMBER_OF_PDFS; ++i) {
447 
448  if (!npe[i].empty() && npe[i].getXmax() > xmax) {
449  xmax = npe[i].getXmax();
450  }
451 
452  }
453 
454  return xmax;
455  }
result_type JFIT::JAbstractRegressor< JShower3EZ , JGandalf >::operator() ( const JShower3EZ value,
T  __begin,
T  __end 
)
inlineinherited

Global fit.

Parameters
valuestart value
__beginbegin of data set
__endend of data set
Returns
chi2

Definition at line 92 of file JRegressor.hh.

93  {
94  static_cast<minimiser_type&>(*this).value = value;
95 
96  return static_cast<minimiser_type&>(*this)(static_cast<regressor_type&>(*this), __begin, __end);
97  }
JModel_t value
Definition: JGandalf.hh:256
JRegressor< JShower3EZ, JGandalf > regressor_type
Definition: JRegressor.hh:79
result_type JFIT::JGandalf< JShower3EZ >::operator() ( const JFunction_t &  fit,
T  __begin,
T  __end,
Args...  args 
)
inlineinherited

Multi-dimensional fit of two data sets.

The fit function should return the equivalent of chi2 for the current value of the model and the given data point as well as the partial derivatives.

Parameters
fitfit function
__beginbegin of data
__endend of data
Returns
chi2

Definition at line 127 of file JGandalf.hh.

128  {
129  using namespace std;
130  using namespace JPP;
131 
132  const int N = parameters.size();
133 
134  H.resize(N);
135  h.resize(N);
136 
137  previous = value;
138 
139  successor.gradient = value;
140  successor.gradient = zero;
141 
142  error = value;
143  error = zero;
144 
145  lambda = LAMBDA_MIN;
146 
147  result_type precessor = result_type(numeric_limits<double>::max(), value);
148 
150 
151  DEBUG("step: " << numberOfIterations << endl);
152 
153  reset();
154 
155  __evaluate__(fit, __begin, __end, args...);
156 
157  DEBUG("lambda: " << SCIENTIFIC(12,5) << lambda << endl);
158  DEBUG("chi2: " << FIXED (12,5) << successor.chi2 << endl);
159 
160  if (successor.chi2 < precessor.chi2) {
161 
162  if (numberOfIterations != 0) {
163 
164  if (fabs(precessor.chi2 - successor.chi2) < EPSILON*fabs(precessor.chi2)) {
165  return successor;
166  }
167 
168  if (lambda > LAMBDA_MIN) {
169  lambda /= LAMBDA_DOWN;
170  }
171  }
172 
173  precessor = successor;
174  previous = value;
175 
176  } else {
177 
178  value = previous;
179  lambda *= LAMBDA_UP;
180 
181  if (lambda > LAMBDA_MAX) {
182  return precessor; // no improvement found
183  }
184 
185  reset();
186 
187  __evaluate__(fit, __begin, __end, args...);
188  }
189 
190  DEBUG("Hesse matrix:" << endl);
191  DEBUG(H << endl);
192 
193 
194  // force definite positiveness
195 
196  for (int i = 0; i != N; ++i) {
197 
198  if (H(i,i) < PIVOT) {
199  H(i,i) = PIVOT;
200  }
201 
202  h[i] = 1.0 / sqrt(H(i,i));
203  }
204 
205 
206  // normalisation
207 
208  for (int i = 0; i != N; ++i) {
209  for (int j = 0; j != i; ++j) {
210  H(j,i) *= h[i] * h[j];
211  H(i,j) = H(j,i);
212  }
213  }
214 
215  for (int i = 0; i != N; ++i) {
216  H(i,i) = 1.0 + lambda;
217  }
218 
219 
220  try {
221  H.invert();
222  }
223  catch (const JException& error) {
224  ERROR("JGandalf: " << error.what() << endl);
225  return precessor;
226  }
227 
228 
229  for (int i = 0; i != N; ++i) {
230 
231  DEBUG("u[" << noshowpos << i << "] = " << showpos << FIXED(15,5) << alias(value, parameters[i]));
232 
233  for (int j = 0; j != N; ++j) {
234  alias(value, parameters[i]) -= H(i,j) * alias(successor.gradient, parameters[j]) * h[i] * h[j];
235  }
236 
237  DEBUG(' ' << FIXED(15,5) << alias(value, parameters[i]) << noshowpos << endl);
238 
239  alias(error, parameters[i]) = h[i];
240  }
241  }
242 
243  return precessor;
244  }
static double EPSILON
maximal distance to minimum
Definition: JGandalf.hh:248
static double PIVOT
minimal value diagonal element of matrix
Definition: JGandalf.hh:253
static const JZero zero
Function object to assign zero value.
Definition: JZero.hh:105
Auxiliary data structure for floating point format specification.
Definition: JPrint.hh:481
std::vector< parameter_type > parameters
Definition: JGandalf.hh:258
void resize(const size_t size)
Resize matrix.
Definition: JMatrixND.hh:359
static double LAMBDA_MIN
minimal value control parameter
Definition: JGandalf.hh:249
static double LAMBDA_DOWN
multiplication factor control parameter
Definition: JGandalf.hh:252
#define ERROR(A)
Definition: JMessage.hh:66
static double alias(const T &model, const typename T::parameter_type parameter)
Read/write access to parameter value by data member.
Definition: JGandalf.hh:326
static double LAMBDA_UP
multiplication factor control parameter
Definition: JGandalf.hh:251
void invert()
Invert matrix according LDU decomposition.
Definition: JMatrixNS.hh:80
void __evaluate__(const JFunction_t &fit, T __begin, T __end, Args...args)
Recursive method for evaluation of fit.
Definition: JGandalf.hh:286
static int MAXIMUM_ITERATIONS
maximal number of iterations
Definition: JGandalf.hh:247
int j
Definition: JPolint.hh:634
static double LAMBDA_MAX
maximal value control parameter
Definition: JGandalf.hh:250
Auxiliary data structure for floating point format specification.
Definition: JPrint.hh:518
then usage $script[input file[working directory[option]]] nWhere option can be N
Definition: JMuonPostfit.sh:37
#define DEBUG(A)
Message macros.
Definition: JMessage.hh:62
std::vector< double > h
Definition: JGandalf.hh:404

Member Data Documentation

JTimeRange JFIT::JRegressor< JShower3EZ, JGandalf >::T_ns
static

Time window with respect to Cherenkov hypothesis [ns].

Definition at line 457 of file JShower3EZRegressor.hh.

double JFIT::JRegressor< JShower3EZ, JGandalf >::Vmax_npe = std::numeric_limits<double>::max()
static

Maximal integral of PDF [npe].

Definition at line 458 of file JShower3EZRegressor.hh.

const int JFIT::JRegressor< JShower3EZ, JGandalf >::NUMBER_OF_PDFS = 2
static

Definition at line 460 of file JShower3EZRegressor.hh.

const JPDFType_t JFIT::JRegressor< JShower3EZ, JGandalf >::pdf_t
static
Initial value:

Definition at line 462 of file JShower3EZRegressor.hh.

PDF.

Definition at line 464 of file JShower3EZRegressor.hh.

M-Estimator function.

Definition at line 466 of file JShower3EZRegressor.hh.

int JFIT::JGandalf< JShower3EZ >::MAXIMUM_ITERATIONS
staticinherited

maximal number of iterations

maximal number of iterations.

Definition at line 247 of file JGandalf.hh.

double JFIT::JGandalf< JShower3EZ >::EPSILON
staticinherited

maximal distance to minimum

maximal distance to minimum.

Definition at line 248 of file JGandalf.hh.

double JFIT::JGandalf< JShower3EZ >::LAMBDA_MIN
staticinherited

minimal value control parameter

Definition at line 249 of file JGandalf.hh.

double JFIT::JGandalf< JShower3EZ >::LAMBDA_MAX
staticinherited

maximal value control parameter

Definition at line 250 of file JGandalf.hh.

double JFIT::JGandalf< JShower3EZ >::LAMBDA_UP
staticinherited

multiplication factor control parameter

Definition at line 251 of file JGandalf.hh.

double JFIT::JGandalf< JShower3EZ >::LAMBDA_DOWN
staticinherited

multiplication factor control parameter

Definition at line 252 of file JGandalf.hh.

double JFIT::JGandalf< JShower3EZ >::PIVOT
staticinherited

minimal value diagonal element of matrix

Definition at line 253 of file JGandalf.hh.

double JFIT::JGandalf< JShower3EZ >::lambda
inherited

Definition at line 255 of file JGandalf.hh.

JShower3EZ JFIT::JGandalf< JShower3EZ >::value
inherited

Definition at line 256 of file JGandalf.hh.

JShower3EZ JFIT::JGandalf< JShower3EZ >::error
inherited

Definition at line 257 of file JGandalf.hh.

Definition at line 258 of file JGandalf.hh.

int JFIT::JGandalf< JShower3EZ >::numberOfIterations
inherited

Definition at line 259 of file JGandalf.hh.

Definition at line 260 of file JGandalf.hh.

template<class T>
int JEEP::JMessage< T >::debug = 0
staticinherited

debug level (default is off).

Definition at line 45 of file JMessage.hh.


The documentation for this struct was generated from the following file: