Rene Brun
Chris Jillings wrote:
>
> Hi Rooters,
> I have found a pathological fit and I am trying to understand why root
> failed to fit properly.
> I have a macro below that starts the fit at a particular set of
> parameter values close to the correct answer. I calcuated the chi square
> two ways. First I looped over the histogram bins and calculated it myself
> in the macro (Answer: 5233). Second, I set very tight bounds on the
> paremeters, effectively forcing migrad not to change them, and then looked
> at the value of FCN from migrad( Answer 5E8!).
>
> The output of the macro is
> (bin #) data model (model-data)**2 total chiSq
> ---------------
> data
>
> :
> :
>
> 70 2138.000000 2790.239235 198.978494 1390.020928
> 71 2547.000000 3419.200307 298.678200 1688.699128
> 72 3267.000000 4196.629537 264.527418 1953.226546
> 73 4010.000000 5158.748900 329.083300 2282.309846
> 74 4942.000000 6350.814741 401.610477 2683.920323
> 75 6087.000000 7829.403125 498.762715 3182.683039
> 76 7695.000000 9665.278362 504.483018 3687.166056
> 77 9704.000000 11946.994351 518.448440 4205.614496
> 78 12155.000000 14785.416850 569.238404 4774.852900
> 79 15641.000000 18319.403236 458.656345 5233.509246
> **********
> ** 1 **SET ERR 1
> **********
> FCN=5.16288e+08 FROM MIGRAD STATUS=CONVERGED 48 CALLS 49
> TOTAL
> EDM=1.25522e-06 STRATEGY= 1 ERROR MATRIX
> ACCURATE
> EXT PARAMETER STEP FIRST
> NO. NAME VALUE ERROR SIZE DERIVATIVE
> 1 1st slope 9.20000e+00 9.23121e-08 2.13212e-02** at limit **
> 2 2nd Amp 2.00000e-01 1.20541e-09 7.70472e-03** at limit **
> 3 2nd slope 5.70000e+00 4.82292e-07 4.87308e-02** at limit **
> 4 scale fac 1.70000e+04 1.64631e-05 9.00421e-03** at limit **
>
>
> So FCN, which from the docs is a chi square, is 5E8, whereas my
> loop calculates 5233.
>
> Any ideas? Thanks once again for your help. The full macro follows.
>
> Chris
>
> {
> gROOT->Reset();
> gStyle->SetOptFit(1);
> gStyle->SetStatX(0.45);
> gStyle->SetStatY(0.85);
>
> TCanvas* c1 = new TCanvas("c1","Ntuple Plots",1); // create canvas
> c1->SetFillColor(kWhite);
>
> TFile* f1 = new TFile("cjj_g10.root","READ");
> TH1F* hi1 = new TH1F("hi1","Angular Resolution: Gamma Energy = 10 MeV
> ",80,-1.0,1.0);
> hi1->SetXTitle("Cos`q#");
> // h509->Draw("Angres>>hi1","Egen>4&&Rfitt<600","goff");
> h503->Draw("(Ue*Uf+Ve*Vf+We*Wf)>>hi1","Rfit<600","goff");
> c1->SetLogy(1);
> TF1* expExp = new TF1("expExp",expExp,-1,1,4);
> expExp->SetParameters(9.2,0.2,5.7,17000);
> expExp->SetParNames("1st slope","2nd Amp","2nd slope","scale fac");
> Double_t x;
> Double_t y;
> Double_t dChiSq;
> Double_t chiSq = 0..0;
> Double_t data;
> Double_t sigSq;
> printf("Hello\n");
> for( Int_t i=0 ; i<80 ; i++) {
> x = -0.9875+i*0.025;
> // printf("%d\t%f\n",i,x);
> y = expExp->Eval(x);
> data = hi1->GetBinContent(i);
> if (data==0.0) sigSq = y;
> else sigSq = data;
> dChiSq = (data-y)*(data-y)/sigSq;
> chiSq += dChiSq;
> printf("%d\t%f\t%f\t%f\t%f\n",i,data,y,dChiSq,chiSq);
> }
> // expExp->Draw();
> expExp->SetParLimits(0,9.2,9.3);
> expExp->SetParLimits(1,0.2,0.21);
> expExp->SetParLimits(2,5.7,5.8);
> expExp->SetParLimits(3,17000.0,17100.0);
>
> TPostScript ps("cjj_g10.eps",114);
> hi1->Fit("expExp","RB");
> // expExp->Draw();
> // hi1->Draw("same");
> c1->Update();
> ps.Close();
>
> }