hello experts,
i am doing a simple fit however is faced with tolerance value. For
func->SetParameter(0,0.00000041)
it confuses MIGRAD. And therefore, if i attempt to retrieve the approximate value from the fit with
GetParameter(0)
it returns me only precision up to 10^-6. That is
GetParameter(0) = 0.000000
Is there a way to rescale this?
The code below…
TH1F *readDataFile(char* fname, int nch) {
int ch,counts;
float fc;
TH1F *h = new TH1F("h1","",nch,-0.5,nch-0.5);
h->Sumw2();
FILE* fp=fopen(fname,"r");
for (int j=0; j<nch; j++) {
fscanf(fp,"%d %f",&ch,&fc);
h->SetBinContent(j+1,fc);
}
fclose(fp);
return h;
}
TH2F* createScan(char* fn, int mi, int ma, int step, int nch) {
char fname[1000];
TH1F *hdum;int z=5;
TH2F *h2 = new TH2F("h2","",nch,-0.5,nch-0.5,((ma-mi)/step)+1,mi-(Double_t)step/2.,ma+(Double_t)step/2.);
for (int i=mi; i<ma+1; i=i+step) {
sprintf(fname,fn,i);
hdum=readDataFile(fname, nch);
for (int j=0; j<nch; j++)
h2->Fill(j,i,hdum->GetBinContent(j+1));
delete hdum;
}
return h2;
}
float pixelAveCount(TH2F *h2, int ichan) {
float pixAveCount;
int nsteps = h2->GetYaxis()->GetNbins();
float pixCount=0.;
for (int istep=1; istep<nsteps+1; istep++) {
pixCount +=h2->GetBinContent (ichan+1, istep);
}
pixAveCount = pixCount/nsteps;
return pixAveCount;
}
bool item_exists(int item, int *array, int arrSize) {
for(int k=0; k < arrSize; k++) {
if(array[k] == item) {
return true;
}
}
return false;
}
void storeToVector(char *fpath0, int chMin, int chMax) {
char fname[256];
char fn[5] = "/run_";
char ext[5] = ".raw";
int begin = 0;
int end = 59;
char specifier[3] = "%d";
float pixelAveCnts;
int filter = 20;
vector<float> avCnt[filter];
int size_goodChan;
Double_t avPixCount;
vector<float> ioc;
char *fbadname ="/home/badChan_AllMods.raw";
char line_bad[80];
int badChan;
vector<int> vec_bad;
int size;
FILE *fbad = fopen(fbadname,"r");
if(fbad==NULL) {
printf("Could not open file %s\n", fbadname);
return;
}
while (fgets(line_bad, 256, fbad) != NULL) {
sscanf(line_bad, "%d",&badChan);
vec_bad.push_back(badChan);
}
size=vec_bad.size();
int *bad=new int[size*sizeof(int)];
for (int ibad =0; ibad<size; ibad++) {
bad[ibad] = vec_bad[ibad];
}
fclose(fbad);
for (int ifilter = 0; ifilter<filter; ifilter++) {
sprintf(fname, "%s%d%s%s%s", fpath0, ifilter,fn, specifier, ext);
TH2F *h2 = createScan(fname, begin,end, 1,128*12*3);
for (int ichan=chMin; ichan<chMax+1; ichan++) {
if(!item_exists(ichan, bad, size)) {
pixelAveCnts = pixelAveCount(h2, ichan);
avCnt[ifilter].push_back(pixelAveCnts);
}else{
;
}
}
delete h2;
}
int igoodChan;
char channel[]= "Good Channel";
char chanNum[100];
vector<float> iocVector = iocData() ;
float iocData_elem;
int ifilter = 19;
Double_t tdead;
Double_t mypar[2], emypar[2];
size_goodChan = avCnt[0].size();
M1:
if(igoodChan>=0 && igoodChan< size_goodChan){
sprintf(chanNum, "%s %d", channel,igoodChan );
TGraph *geff = new TGraph();
geff ->SetLineColor(2);
geff ->SetLineWidth(3);
geff ->SetTitle(chanNum);
for (int ifilter =0; ifilter<filter; ifilter++) {
avPixCount = avCnt[ifilter][igoodChan];
iocData_elem = iocVector[ifilter];
geff->SetPoint(ifilter, iocData_elem, avPixCount);
}
tdead = fitParalyzable(geff,mypar, emypar, 1);
geff->Draw("AL*");
c1->Modified();
c1->Update();
igoodChan++;
if (igoodChan<=size_goodChan)
goto M1;
}
return;
}
Double_t paralyzable(Double_t *x, Double_t *par) {
Double_t rate;
rate = x[0]*par[1]*TMath::Exp((-1)*x[0]*par[0]*par[1]);
return rate;
}
Double_t fitParalyzable(TGraph *geff, Double_t *mypar, Double_t *emypar, int plot) {
Double_t tdead;
Double_t efficiency;
float chi;
Double_t min=0.;
Double_t max=3900000.;
TF1 *func = new TF1("func",paralyzable,min, max,2);
func->SetParameters(mypar);
func->SetParameter(0,0.0000004); // i am faced with troubles here. For Double_t for the second argument, how does one re-scale this? I am limited to 10^-6 for Double_t if I retrieve back the approximate value from the fit.....and consequently it screws the MIGRAD
func->SetParLimits(0,0.,1.);
func->SetParameter(1,0.91);
func->SetParLimits(0,0.,1);
func->SetParNames("par0", "par1");
func->SetRange(min,max);
if (plot) {
geff->Fit("func","R");
}else{
geff->Fit("func","R0QV");
}
TF1* fitfunc;
fitfunc = geff->GetFunction("func");
fitfunc->GetParameters(mypar);
emypar[0]=fitfunc->GetParError(0);
emypar[1]=fitfunc->GetParError(1);
// chi = fitfunc->GetChisquare()/(fitfunc->GetNDF());
tdead = fitfunc->GetParameter(0);
efficiency = fitfunc->GetParameter(1);
delete func;
if (plot)
;
return tdead;
}