forked from amissert/atmFitTools
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathhistoCompare.h
144 lines (130 loc) · 4.93 KB
/
histoCompare.h
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
#ifndef HISTCOMPARE_H
#define HISTCOMPARE_H
#include "shared.h"
#include "histoManager.h"
#include "TMath.h"
#include "TRandom2.h"
#include "histoTransforms.cxx"
#include "TFitter.h"
#include "TCanvas.h"
#include "TGraphAsymmErrors.h"
#include "likelihood.cxx"
using namespace std;
//class to compare histograms and evaluate likelihoods
class histoCompare{
public:
//constructors//
histoCompare(const char* parfile, bool sepearteneutmode = false); //construct from parameter file
histoCompare(); //standard constructor
//internal variables
bool separateNeutMode;
TString nameTag; //name associated with this instance
int nSamp; //number of samples
int nBin; //number of bins
int nComp; //number of components
int nAtt; //nummboer of attributes
int nMode;
double tunePar; //tuning parameter for MCMC
//tools for histogram manager management
//created histo manager from file
void readFromFile(const char* rootname,int isamp,int ibin, int icomp, int natt);
void readFromFile(const char* rootname,int isamp,int ibin, int icomp, int imode, int natt);
histoManager* hManager;
//atmospheric pars
atmFitPars* thePars;
sharedPars* runPars;
int MCMCNSteps;
double Norm;
double Par[NBINMAX][NCOMPMAX][NATTMAX][2];
double sysPar[NSYSPARMAX];
double sysParUnc[NSYSPARMAX];
double fixPar[NBINMAX][NCOMPMAX][NATTMAX][2];
double bestPar[NBINMAX][NCOMPMAX][NATTMAX][2];
// double errParLo[NBINMAX][NCOMPMAX][NATTMAX][2];
// double errParHi[NBINMAX][NCOMPMAX][NATTMAX][2];
double errParLo[1000];
double errParHi[1000];
TString parName[NBINMAX][NCOMPMAX][NATTMAX][2];
TString binName[NBINMAX];
TString compName[NCOMPMAX];
TString attName[NCOMPMAX];
void setBinName(int ibin, const char* name){binName[ibin]=name;}
void setCompName(int icomp, const char* name){compName[icomp]=name;}
void setAttName(int iatt, const char* name){attName[iatt]=name;}
void setupPars(int nsyspars=0); //sets up all parameters
void setupPars(atmFitPars *a);
//post-fit toolts
void profileL(int ibin, int icomp, int iatt, int imod, double range, int npts=100);
void profileL(int ipar,double range, int npts=100,int sameflg=0);
void showFitHisto(int isamp,int ibin,int icomp,int iatt);
void showFitEffect(int isamp,int ibin,int icomp,int iatt);
void showFitResult(int isamp,int ibin,int iatt);
void showFitPars(int ibin,int iatt,int imod);
void showModHiso(int isamp,int ibin, int icomp, int iatt, double smear, double bias);
void runMCMC(int nsteps);
// double getErrLo(int ibin,int icomp,int iatt,int imod);
double getErrLo(int isyst);
double getErrHi(int isyst);
TH1D* getModifiedHisto(int ibin, int icomp, int iatt){return hManager->getSumHistogramMod(ibin,icomp,iatt);}
// TH1D* hMod[NSAMPMAX][NBINMAX][NCOMPMAX][NATTMAX];
//initialize all necessary components
void initialize(histoManager* hm, atmFitPars* apars);
void setupSplines(const char* fname,int npar){hManager->readSplinesFromFile(fname,npar);}
//tools for adding histogram directly..for debugging
void addHistogram(TH1D* h,int dataflg);
int rebinFactor;
void setRebinFactor(int ifact){rebinFactor=ifact;}
TH1D* hData[10];
TH1D* hMC[10];
TH1D* hModDebug;
TH1D* hMod;
TH1D* hTmp; //temporary histogram container
TH1D* hPar; //for fit parameters
TGraphAsymmErrors* gPar;
TH1D* hParErrLo;
TH1D* hParErrHi;
TH1D* hProf[2];
TH1D* showSmear(TH1D* h, double smear, double bias);
void showMod(int imchist);
TH1D* hTot;
TCanvas* cc;
int nDataHist;
int nMCHist;
int useLnLType;
double parDebug[10][2];
// double getLDebug();
double cScale; //correction scale for likelihood
//////////////////////////////////////////////////////////
//flags
int flgFixAllSmearPars;
//likelihood evaluateions
double getSumSq(TH1D* h1, TH1D* h2);
double getLnL(TH1D* h1, TH1D* h2,double hnorm = 1.);
double getNDiff();
double getTotSumSq();
double getTotLnL();
static void sumSqWrapper(int& ndim, double* gout, double& result, double par[], int flg);
static void lnLWrapper(int& ndim, double* gout, double& result, double par[], int flg);
void getTotLnL1D(double& result,int npar, double par[]);
//for debuggint and play
double getTotSumSqDebug();
static void sumSqDebugWrapper(int& ndim, double* gout, double& result, double par[], int flg);
static void nDiffDebugWrapper(int& ndim, double* gout, double& result, double par[], int flg);
void minSumSqDebug();
void minSumSq();
void sumSqPrefit();
void LnLFit();
void LnLPreFit();
void singleParFit(int ipar);
void sysParFit();
void drawResult(int ihist);
void timetest(int ntry);
void saveFitPars(const char* filename); //< write parameters and outputs to a file
void readFitPars(const char* filename); //< read parameters from a file
void tuneMCMC(int ncyles=1,int nsteps=150,double goal=0.25);
void tuneMCMCOld(int ncyles=1,int nsteps=150,double goal=0.25);
void calcRoughParErr();
//staticthis for fits
static histoCompare* staticthis;
};
#endif