28 using TMath::PoissonI;
39 inline double g1(
const JGauss&
gauss,
const JElement_t& point)
41 const double u = (point.getX() - gauss.
mean) / gauss.
sigma;
42 const double fs = gauss.
signal *
exp(-0.5*u*u) / (sqrt(2.0*
PI) * gauss.
sigma);
44 const double f1 = fs + fb;
46 return -
log(PoissonI(point.getY(),
f1));
57 int main(
int argc,
char **argv)
70 JParser<> zap(
"Program to test JSimplex algorithm.");
80 catch(
const exception& error) {
81 FATAL(error.what() << endl);
84 ASSERT(numberOfEvents > 0);
88 TF1 fs(
"fs",
"exp(-0.5 * (x-[0])*(x-[0]) / ([1]*[1]))");
91 fs.FixParameter(0, gauss.
mean);
92 fs.FixParameter(1, gauss.
sigma);
95 const Double_t
xmin = -5.0;
96 const Double_t
xmax = +5.0;
103 TH1D
H[] = { TH1D(
"ha",
"", 101, -0.1, +0.1),
104 TH1D(
"hb",
"", 101, -0.1, +0.1),
105 TH1D(
"hc",
"", 101, -100.0, +100.0),
106 TH1D(
"hd",
"", 101, -100.0, +100.0) };
110 for (
int i = 0; i != numberOfEvents; ++i) {
112 STATUS(
"event: " << setw(10) << i <<
'\r');
DEBUG(endl);
114 TH1D h0(
"h0", NULL, nx,
xmin,
xmax);
116 h0.FillRandom(
"fs", (Int_t) gauss.
signal);
117 h0.FillRandom(
"fb", (Int_t) gauss.
background);
121 for (Int_t i = 1; i <= h0.GetNbinsX(); ++i) {
122 data.push_back(JElement_t(h0.GetBinCenter (i),
123 h0.GetBinContent(i)));
137 h0.GetEntries() * (
xmax -
xmin) / nx - h0.GetMinimum(),
142 for (
size_t i = 0; i != fit.
step.size(); ++i) {
143 DEBUG(
"Step size " << i <<
' ' << fit.
step[i] << endl);
148 const double chi2 = fit(
g1, data.begin(), data.end());
153 DEBUG(
"Chi2 " << chi2 << endl);
160 for (
int i = 0; i !=
sizeof(
Q)/
sizeof(
Q[0]); ++i) {
166 for (
int i = 0; i !=
sizeof(
Q)/
sizeof(
Q[0]); ++i) {
171 timer.print(cout,
true,
micro_t);
178 for (
int i = 0; i !=
sizeof(
H)/
sizeof(
H[0]); ++i) {
186 for (
int i = 0; i !=
sizeof(
Q)/
sizeof(
Q[0]); ++i) {
190 ASSERT(
Q[0].getSTDev() < precision.mean);
191 ASSERT(
Q[1].getSTDev() < precision.sigma);
192 ASSERT(
Q[2].getSTDev() < precision.signal);
193 ASSERT(
Q[3].getSTDev() < precision.background);
Utility class to parse command line options.
Q(UTCMax_s-UTCMin_s)-livetime_s
int main(int argc, char *argv[])
The elements in a collection are sorted according to their abscissa values and a given distance opera...
static const double H
Planck constant [eV s].
then fatal Wrong number of arguments fi set_variable STRING $argv[1] set_variable DETECTORXY_TXT $WORKDIR $DETECTORXY_TXT tail read X Y CHI2 RMS printf optimum n $X $Y $CHI2 $RMS awk v Y
double getMean(vector< double > &v)
get mean of vector content
const JPolynome f1(1.0, 2.0, 3.0)
Function.
#define ASSERT(A,...)
Assert macro.
I/O formatting auxiliaries.
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Auxiliary class for CPU timing and usage.
static const double PI
Mathematical constants.
std::ostream & shortprint(std::ostream &out)
Set short printing.
General purpose messaging.
std::ostream & longprint(std::ostream &out)
Set long printing.
Simple fit method based on Powell's algorithm, see reference: Numerical Recipes in C++...
Utility class to parse command line options.
then set_variable NUMBER_OF_TESTS else set_variable NUMBER_OF_TESTS fi function gauss()
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
then cat $TRIPOD_INITIAL<< EOF1 256877.5 4743716.7-2438.42 256815.5 4743395.0-2435.53 257096.2 4743636.0-2439.5EOFfiif[[!-f $DETECTOR]];then JEditDetector-a $DETECTOR_INITIAL-s"-1 addz -6.9"-o $DETECTOR--!eval`JPrintDetector-a $DETECTOR-O SUMMARY`for STRING in ${STRINGS[*]};do set_variable MODULE`getModule-a $DETECTOR-L"$STRING 0"`JEditDetector-a $DETECTOR-M"$MODULE setz -2.9"-o $DETECTOR--!donefiif[[!-f $TRIPOD]];then cp-p $TRIPOD_INITIAL $TRIPODfiJAcoustics.sh $DETECTOR_IDcat > acoustics_trigger_parameters txt<< EOFQ=0.0;TMax_s=0.020;numberOfHits=90;EOFJAcousticsEventBuilder.sh $DETECTOR $RUNS[*]INPUT_FILES=(`ls KM3NeT_ ${(l:8::0::0:) DETECTOR_ID}_0 *${^RUNS}_event.root`) cd $WORKDIRif[!$HOMEDIR-ef $WORKDIR];then cp-p $HOMEDIR/$DETECTOR $WORKDIR cp-p $HOMEDIR/${^ACOUSTICS_KEYS}.txt $WORKDIR cp-p $HOMEDIR/${^INPUT_FILES}$WORKDIRfisource $JPP_DIR/examples/JAcoustics/acoustics-fit-toolkit.shtimer_startinitialise stage_b 1 0 100.0e-6 0.002 0.1 0 > &stage log
then fatal Wrong number of arguments fi set_variable DETECTOR $argv[1] set_variable STRING $argv[2] set_array QUANTILES set_variable FORMULA *[0] exp(-0.5 *(x-[1])*(x-[1])/([2]*[2]))" set_variable MODULE `getModule -a $DETECTOR -L "$STRING 0"` source JAcousticsToolkit.sh typeset -A TRIPODS get_tripods $WORKDIR/tripod.txt TRIPODS XMEAN
#define DEBUG(A)
Message macros.
std::vector< JModel_t > step
Double_t g1(const Double_t x)
Function.