1 #ifndef __JFIT__JGRADIENT__
2 #define __JFIT__JGRADIENT__
21 namespace JPP {
using namespace JFIT; }
40 virtual void apply(
const double step) = 0;
48 public std::shared_ptr<JParameter_t>
91 const int debug = 3) :
123 return numeric_limits<double>::max();
128 const size_t N = this->size();
133 for (
size_t i = 0;
i !=
N; ++
i) {
152 for (
double ds = 1.0; ds > 1.0e-3; ) {
154 this->
move(+1.0 * ds);
158 DEBUG(
"chi2[3] " << setw(4) << m <<
' ' <<
FIXED(12,5) << chi2[3] <<
' ' <<
FIXED(12,5) << ds << endl);
160 if (chi2[3] < chi2[2]) {
184 for ( ; m != 0; --m) {
185 this->
move(-1.0 * ds);
192 this->
move(-1.0 * ds);
194 if (chi2[2] < chi2[3]) {
202 const double f21 = chi2[2] - chi2[1];
203 const double f23 = chi2[2] - chi2[3];
205 const double xs = 0.5 * (f21 - f23) / (f23 + f21);
207 this->
move(+1.0 * xs * ds);
211 if (chi2[3] < chi2[2]) {
217 this->
move(-1.0 * xs * ds);
232 if (fabs(chi2[2] - chi2[0]) <
epsilon * 0.5 * (fabs(chi2[0]) + fabs(chi2[2]))) {
244 for (
size_t i = 0;
i !=
N; ++
i){
255 for (
size_t i = 0;
i !=
N; ++
i){
257 H[
i] = G[
i] + dgg * H[
i];
287 const size_t N = this->size();
299 for (
size_t i = 0;
i !=
N; ++
i) {
300 if ((*
this)[
i].
name.size() > width) {
301 width = (*this)[
i].name.size();
307 for (
size_t i = 0;
i !=
N; ++
i) {
309 if ((*
this)[
i].value != 0.0) {
311 (*this)[
i]->apply(+0.5 * (*
this)[
i].value);
315 (*this)[
i]->apply(-0.5 * (*
this)[
i].value);
316 (*this)[
i]->apply(-0.5 * (*
this)[
i].value);
322 (*this)[
i]->apply(+0.5 * (*
this)[
i].value);
334 DEBUG(setw(width) << left <<
"|gradient|" << right <<
' ' <<
FIXED(12,5) << sqrt(
V) << endl);
348 for (
size_t i = 0;
i != this->size(); ++
i) {
349 (*this)[
i ]->apply((*
this)[
i ].value *
gradient[
i ] * factor);
351 }
else if (factor < 0.0) {
352 for (
size_t i = this->size();
i != 0; --
i) {
353 (*this)[
i-1]->apply((*
this)[
i-1].value *
gradient[
i-1] * factor);
JModifier_t(const std::string &name, JParameter_t *parameter, const double value)
Constructor.
void move(const double factor)
Move.
virtual ~JParameter_t()
Virtual destructor.
then echo Enter input within $TIMEOUT_S seconds echo n User name
double operator()(const T &getChi2)
Fit.
static const double H
Planck constant [eV s].
Auxiliary data structure for floating point format specification.
V(JDAQEvent-JTriggerReprocessor)*1.0/(JDAQEvent+1.0e-10)
size_t Nextra
maximum number of extra steps
do set_variable OUTPUT_DIRECTORY $WORKDIR T
JGradient(const size_t Nmax=std::numeric_limits< size_t >::max(), const size_t Nextra=0, const double epsilon=1.0e-4, const int debug=3)
Constructor.
std::vector< double > gradient
virtual void apply(const double step)=0
Apply step.
size_t Nmax
maximum number of iterations
Auxiliary data structure for fit parameter.
General purpose messaging.
Auxiliary data structure for editable parameter.
double evaluate(const T &getChi2)
Evaluate gradient.
then usage $script< input file >[option[primary[working directory]]] nWhere option can be N
then if[[!-f $DETECTOR]] then JDetector sh $DETECTOR fi cat $WORKDIR trigger_parameters txt<< EOFtrigger3DMuon.enabled=1;trigger3DMuon.numberOfHits=5;trigger3DMuon.gridAngle_deg=1;ctMin=0.0;TMaxLocal_ns=15.0;EOF set_variable TRIGGEREFFICIENCY_TRIGGERED_EVENTS_ONLY INPUT_FILES=() for((i=1;$i<=$NUMBER_OF_RUNS;++i));do JSirene.sh $DETECTOR $JPP_DATA/genhen.km3net_wpd_V2_0.evt.gz $WORKDIR/sirene_ ${i}.root JTriggerEfficiency.sh $DETECTOR $DETECTOR $WORKDIR/sirene_ ${i}.root $WORKDIR/trigger_efficiency_ ${i}.root $WORKDIR/trigger_parameters.txt $JPP_DATA/PMT_parameters.txt INPUT_FILES+=($WORKDIR/trigger_efficiency_ ${i}.root) done for ANGLE_DEG in $ANGLES_DEG[*];do set_variable SIGMA_NS 3.0 set_variable OUTLIERS 3 set_variable OUTPUT_FILE $WORKDIR/matrix\[${ANGLE_DEG}\deg\].root $JPP_DIR/examples/JReconstruction-f"$INPUT_FILES[*]"-o $OUTPUT_FILE-S ${SIGMA_NS}-A ${ANGLE_DEG}-O ${OUTLIERS}-d ${DEBUG}--!fiif[[$OPTION=="plot"]];then if((0));then for H1 in h0 h1;do JPlot1D-f"$WORKDIR/matrix["${^ANGLES_DEG}" deg].root:${H1}"-y"1 2e3"-Y-L TR-T""-\^"number of events [a.u.]"-> o chi2
size_t numberOfIterations
double getChi2(const double P)
Get chi2 corresponding to given probability.
Auxiliary data structure for floating point format specification.
#define DEBUG(A)
Message macros.