File indexing completed on 2024-04-06 11:59:32
0001 #define analyse_residuals_cxx
0002 #include "analyse_residuals.h"
0003 #include <TH1F.h>
0004 #include <TH2.h>
0005 #include <TStyle.h>
0006 #include <TCanvas.h>
0007 #include <iostream>
0008 #include <TF1.h>
0009 #include <TMath.h>
0010 #include "TFile.h"
0011
0012 double ng(Double_t* xx, Double_t* par)
0013 {
0014 Double_t x = xx[0];
0015 Double_t norm = par[0];
0016 Double_t mean = par[1];
0017 Double_t sigma = par[2];
0018
0019 return norm * 1.0/sqrt(2.0*TMath::Pi())/sigma * exp(-0.5*(x-mean)*(x-mean)/sigma/sigma);
0020 }
0021
0022 const double a_min = 1.37078;
0023 const double a_max = 1.77078;
0024 const double a_bin = 0.10000;
0025
0026 double aa_min[3] = {1.50, 1.45, 1.40};
0027 double aa_max[3] = {1.75, 1.70, 1.65};
0028
0029 double ys_bl[6];
0030 double ys_bh[6];
0031
0032 TCanvas* can_y_barrel_sizey_alpha[6][4];
0033 TCanvas* can_y_barrel_sizey[6];
0034 TH1F* h_yres_npix_alpha_beta[6][4][10];
0035 TH1F* h_yres_npix_alpha[6][4];
0036 TH1F* h_yres_npix_alpha_rms[6][4];
0037
0038 TCanvas* can_x_barrel_sizex_beta[3][4];
0039 TCanvas* can_x_barrel_sizex[3];
0040 TH1F* h_xres_npix_beta_alpha[3][4][10];
0041 TH1F* h_xres_npix_beta[3][4];
0042 TH1F* h_xres_npix_beta_rms[3][4];
0043
0044
0045
0046
0047
0048
0049
0050
0051
0052
0053
0054
0055
0056
0057
0058
0059
0060 TCanvas* can_y_forward_sizey[2];
0061 TCanvas* can_y_forward;
0062 TH1F* h_forward_yres_npix_beta[2][10];
0063 TH1F* h_forward_yres_npix[2];
0064 TH1F* h_forward_yres_npix_rms[2];
0065
0066 TCanvas* can_x_forward_sizex[2];
0067 TCanvas* can_x_forward;
0068 TH1F* h_forward_xres_npix_alpha[2][10];
0069 TH1F* h_forward_xres_npix[2];
0070 TH1F* h_forward_xres_npix_rms[2];
0071
0072 const char* fname;
0073
0074 void analyse_residuals::Loop()
0075 {
0076 if (fChain == 0) return;
0077
0078 Long64_t nentries = fChain->GetEntries();
0079 cout << "nentries = " << nentries << endl;
0080
0081 bool do_residuals = true;
0082 bool do_plots = true;
0083
0084 bool do_yb = true;
0085 bool do_xb = true;
0086 bool do_yf = true;
0087 bool do_xf = true;
0088
0089 bool do_fix = true;
0090
0091 if ( do_residuals )
0092 fname = "residuals.dat";
0093 else
0094 fname = "pulls.dat";
0095
0096 printf("%s \n", fname);
0097 FILE* datfile;
0098 if ( (datfile = fopen(fname, "w"))==NULL )
0099 { cout << "could not open the output file." << endl;
0100 exit(-1);
0101 }
0102
0103 char hname[100];
0104
0105 ys_bl[0] = 0.05;
0106 ys_bh[0] = 0.50;
0107
0108 ys_bl[1] = 0.15;
0109 ys_bh[1] = 0.90;
0110
0111 ys_bl[2] = 0.70;
0112 ys_bh[2] = 1.05;
0113
0114 ys_bl[3] = 0.95;
0115 ys_bh[3] = 1.15;
0116
0117 ys_bl[4] = 1.15;
0118 ys_bh[4] = 1.20;
0119
0120 ys_bl[5] = 1.20;
0121 ys_bh[5] = 1.40;
0122
0123 if ( do_yb )
0124 {
0125 for (int i=0; i<6; ++i)
0126 {
0127 for (int j=0; j<4; ++j)
0128 {
0129 sprintf(hname, "can_y_barrel_sizey_alpha_%i_%i", i, j);
0130 can_y_barrel_sizey_alpha[i][j] = new TCanvas(hname, hname, 1200, 500);
0131 can_y_barrel_sizey_alpha[i][j]->Divide(5,2);
0132 }
0133
0134 sprintf(hname, "can_y_barrel_sizey_alpha_%i", i);
0135 can_y_barrel_sizey[i] = new TCanvas(hname, hname, 1200, 350);
0136 can_y_barrel_sizey[i]->Divide(4);
0137 }
0138 }
0139
0140 if ( do_xb )
0141 {
0142 for (int i=0; i<3; ++i)
0143 {
0144 for (int j=0; j<4; ++j)
0145 {
0146 sprintf(hname, "can_x_barrel_sizex_beta_%i_%i", i, j);
0147 can_x_barrel_sizex_beta[i][j] = new TCanvas(hname, hname, 1200, 500);
0148 can_x_barrel_sizex_beta[i][j]->Divide(5,2);
0149 }
0150
0151 sprintf(hname, "can_x_barrel_sizex_beta_%i", i);
0152 can_x_barrel_sizex[i] = new TCanvas(hname, hname, 1200, 350);
0153 can_x_barrel_sizex[i]->Divide(4);
0154 }
0155 }
0156
0157 if ( do_yf )
0158 {
0159 for (int i=0; i<2; ++i)
0160 {
0161 sprintf(hname, "can_y_forward_sizey_%i", i);
0162 can_y_forward_sizey[i] = new TCanvas(hname, hname, 1200, 500);
0163 can_y_forward_sizey[i]->Divide(5,2);
0164 }
0165
0166 sprintf(hname, "can_y_forward");
0167 can_y_forward = new TCanvas(hname, hname, 800, 400);
0168 can_y_forward->Divide(2);
0169 }
0170
0171 if ( do_xf )
0172 {
0173 for (int i=0; i<2; ++i)
0174 {
0175 sprintf(hname, "can_x_forward_sizex_%i", i);
0176 can_x_forward_sizex[i] = new TCanvas(hname, hname, 1200, 500);
0177 can_x_forward_sizex[i]->Divide(5,2);
0178 }
0179
0180 sprintf(hname, "can_x_forward");
0181 can_x_forward = new TCanvas(hname, hname, 800, 400);
0182 can_x_forward->Divide(2);
0183 }
0184
0185
0186 for (int i=0; i<6; ++i)
0187 for (int j=0; j<4; ++j)
0188 for (int k=0; k<10; ++k)
0189 {
0190 sprintf(hname, "h_yres_npix_alpha_beta_%i_%i_%i", i, j, k );
0191 if ( do_residuals )
0192 h_yres_npix_alpha_beta[i][j][k]= new TH1F(hname, hname, 100, -0.02, 0.02);
0193 else
0194 h_yres_npix_alpha_beta[i][j][k]= new TH1F(hname, hname, 100, -10.0, 10.0);
0195 }
0196
0197 for (int i=0; i<6; ++i)
0198 for (int j=0; j<4; ++j)
0199 {
0200 sprintf(hname, "h_yres_npix_alpha_%i_%i", i, j );
0201 h_yres_npix_alpha[i][j] = new TH1F(hname, hname, 10, ys_bl[i], ys_bh[i]);
0202
0203 sprintf(hname, "h_yres_npix_alpha_rms_%i_%i", i, j );
0204 h_yres_npix_alpha_rms[i][j] = new TH1F(hname, hname, 10, ys_bl[i], ys_bh[i]);
0205 h_yres_npix_alpha_rms[i][j]->SetLineColor(kRed);
0206 if ( do_residuals )
0207 {
0208 h_yres_npix_alpha_rms[i][j]->SetMinimum(0.0);
0209 h_yres_npix_alpha_rms[i][j]->SetMaximum(0.0060);
0210 }
0211 else
0212 {
0213 h_yres_npix_alpha_rms[i][j]->SetMinimum(0.0);
0214 h_yres_npix_alpha_rms[i][j]->SetMaximum(2.0);
0215 }
0216 }
0217
0218 for (int i=0; i<3; ++i)
0219 for (int j=0; j<4; ++j)
0220 for (int k=0; k<10; ++k)
0221 {
0222 sprintf(hname, "h_xres_npix_beta_alpha_%i_%i_%i", i, j, k );
0223 if ( do_residuals )
0224 h_xres_npix_beta_alpha[i][j][k] = new TH1F(hname, hname, 100, -0.01, 0.01);
0225 else
0226 h_xres_npix_beta_alpha[i][j][k] = new TH1F(hname, hname, 100, -10.0, 10.0);
0227
0228
0229
0230
0231
0232
0233
0234
0235
0236
0237
0238
0239
0240
0241
0242 }
0243
0244 for (int i=0; i<3; ++i)
0245 for (int j=0; j<4; ++j)
0246 {
0247 sprintf(hname, "h_xres_npix_beta_%i_%i", i, j );
0248 h_xres_npix_beta[i][j]= new TH1F(hname, hname, 10, a_min, a_max);
0249
0250 sprintf(hname, "h_xres_npix_beta_rms_%i_%i", i, j );
0251 h_xres_npix_beta_rms[i][j]= new TH1F(hname, hname, 10, a_min, a_max);
0252 h_xres_npix_beta_rms[i][j]->SetLineColor(kRed);
0253 if ( do_residuals )
0254 {
0255 h_xres_npix_beta_rms[i][j]->SetMinimum(0.0);
0256 h_xres_npix_beta_rms[i][j]->SetMaximum(0.0060);
0257 }
0258 else
0259 {
0260 h_xres_npix_beta_rms[i][j]->SetMinimum(0.0);
0261 h_xres_npix_beta_rms[i][j]->SetMaximum(2.0);
0262 }
0263
0264
0265
0266
0267
0268
0269
0270
0271
0272
0273
0274
0275
0276
0277
0278
0279
0280
0281
0282
0283
0284
0285
0286
0287
0288
0289
0290
0291
0292
0293
0294
0295
0296
0297
0298
0299
0300 }
0301
0302
0303 for (int i=0; i<2; ++i)
0304 for (int j=0; j<10; ++j)
0305 {
0306 sprintf(hname, "h_forward_yres_npix_beta_%i_%i", i, j );
0307 if ( do_residuals )
0308 h_forward_yres_npix_beta[i][j] = new TH1F(hname, hname, 100, -0.01, 0.01);
0309 else
0310 h_forward_yres_npix_beta[i][j] = new TH1F(hname, hname, 100, -10.0, 10.0);
0311 }
0312
0313 for (int i=0; i<2; ++i)
0314 {
0315 sprintf(hname, "h_forward_yres_npix_%i", i );
0316 h_forward_yres_npix[i] = new TH1F(hname, hname, 10, 0.3, 0.4);
0317
0318 sprintf(hname, "h_forward_yres_npix_rms_%i", i );
0319 h_forward_yres_npix_rms[i] = new TH1F(hname, hname, 10, 0.3, 0.4);
0320 h_forward_yres_npix_rms[i]->SetLineColor(kRed);
0321 if ( do_residuals )
0322 {
0323 h_forward_yres_npix_rms[i]->SetMinimum(0.0);
0324 h_forward_yres_npix_rms[i]->SetMaximum(0.0040);
0325 }
0326 else
0327 {
0328 h_forward_yres_npix_rms[i]->SetMinimum(0.0);
0329 h_forward_yres_npix_rms[i]->SetMaximum(2.0);
0330 }
0331 }
0332
0333 for (int i=0; i<2; ++i)
0334 for (int j=0; j<10; ++j)
0335 {
0336 sprintf(hname, "h_forward_xres_npix_alpha_%i_%i", i, j );
0337 if ( do_residuals )
0338 h_forward_xres_npix_alpha[i][j] = new TH1F(hname, hname, 100, -0.01, 0.01);
0339 else
0340 h_forward_xres_npix_alpha[i][j] = new TH1F(hname, hname, 100, -10.0, 10.0);
0341 }
0342
0343 sprintf(hname, "h_forward_xres_npix_%i", 0 );
0344 h_forward_xres_npix[0] = new TH1F(hname, hname, 10, 0.15, 0.3);
0345 sprintf(hname, "h_forward_xres_npix_%i", 1 );
0346 h_forward_xres_npix[1] = new TH1F(hname, hname, 10, 0.15, 0.5);
0347
0348 sprintf(hname, "h_forward_xres_npix_rms_%i", 0 );
0349 h_forward_xres_npix_rms[0] = new TH1F(hname, hname, 10, 0.15, 0.3);
0350 h_forward_xres_npix_rms[0]->SetLineColor(kRed);
0351 sprintf(hname, "h_forward_xres_npix_rms_%i", 1 );
0352 h_forward_xres_npix_rms[1] = new TH1F(hname, hname, 10, 0.15, 0.5);
0353 h_forward_xres_npix_rms[1]->SetLineColor(kRed);
0354
0355 if ( do_residuals )
0356 {
0357 h_forward_xres_npix_rms[0]->SetMinimum(0.0);
0358 h_forward_xres_npix_rms[0]->SetMaximum(0.0040);
0359 h_forward_xres_npix_rms[1]->SetMinimum(0.0);
0360 h_forward_xres_npix_rms[1]->SetMaximum(0.0040);
0361 }
0362 else
0363 {
0364 h_forward_xres_npix_rms[0]->SetMinimum(0.0);
0365 h_forward_xres_npix_rms[0]->SetMaximum(2.0);
0366 h_forward_xres_npix_rms[1]->SetMinimum(0.0);
0367 h_forward_xres_npix_rms[1]->SetMaximum(2.0);
0368 }
0369
0370 for (Long64_t jentry=0; jentry<nentries; jentry++)
0371 {
0372 Long64_t ientry = LoadTree(jentry);
0373 if (ientry < 0) break;
0374 fChain->GetEntry(jentry);
0375
0376 if ( !do_residuals )
0377 {
0378 rechitresx = rechitpullx;
0379 rechitresy = rechitpully;
0380 }
0381
0382 if ( do_fix )
0383 {
0384 alpha = trk_alpha;
0385 beta = trk_beta ;
0386 }
0387
0388 double alpha_rad = fabs(alpha)/180.0*TMath::Pi();
0389 double beta_rad = fabs(beta) /180.0*TMath::Pi();
0390 double betap_rad = fabs( TMath::Pi()/2.0 - beta_rad );
0391 double alphap_rad = fabs( TMath::Pi()/2.0 - alpha_rad );
0392
0393 if ( subdetId == 1 )
0394 {
0395
0396 int sizey = nypix;
0397
0398
0399 if ( !( sizey == 1 && bigy == 1 ) )
0400 {
0401 if ( sizey > 6 ) sizey = 6;
0402
0403 int ind_sizey = sizey - 1;
0404 int ind_alpha = -9999;
0405 int ind_beta = -9999;
0406
0407 if ( alpha_rad <= a_min ) ind_alpha = 0;
0408 else if ( alpha_rad >= a_max ) ind_alpha = 3;
0409 else if ( alpha_rad > a_min &&
0410 alpha_rad < a_max )
0411 {
0412 double binw = ( a_max - a_min ) / 4.0;
0413 ind_alpha = (int)( ( alpha_rad - a_min ) / binw );
0414 }
0415 else cout << " Wrong alpha_rad = " << alpha_rad << endl << endl;
0416
0417 if ( betap_rad <= ys_bl[sizey-1] ) ind_beta = 0;
0418 else if ( betap_rad >= ys_bh[sizey-1] ) ind_beta = 9;
0419 else if ( betap_rad > ys_bl[sizey-1] &&
0420 betap_rad < ys_bh[sizey-1] )
0421 {
0422 double binw = ( ys_bh[sizey-1] - ys_bl[sizey-1] ) / 8.0;
0423 ind_beta = 1 + (int)( ( betap_rad - ys_bl[sizey-1] ) / binw );
0424 }
0425 else cout << " Wrong betap_rad = " << betap_rad << endl << endl;
0426
0427 h_yres_npix_alpha_beta[ind_sizey][ind_alpha][ind_beta]->Fill( rechitresy );
0428
0429 }
0430
0431
0432 int sizex = nxpix;
0433
0434 if ( !( sizex == 1 && bigx == 1 ) )
0435 {
0436 if ( sizex > 3 ) sizex = 3;
0437
0438 int ind_sizex = sizex - 1;
0439 int ind_beta = -9999;
0440 int ind_alpha = -9999;
0441
0442 if ( betap_rad <= 0.7 ) ind_beta = 0;
0443 else if ( 0.7 < betap_rad && betap_rad <= 1.0 ) ind_beta = 1;
0444 else if ( 1.0 < betap_rad && betap_rad <= 1.2 ) ind_beta = 2;
0445 else if ( 1.2 <= betap_rad ) ind_beta = 3;
0446 else cout << " Wrong betap_rad = " << betap_rad << " !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! " << endl;
0447
0448 if ( alpha_rad <= aa_min[ind_sizex] ) ind_alpha = 0;
0449 else if ( alpha_rad >= aa_max[ind_sizex] ) ind_alpha = 9;
0450 else
0451 ind_alpha = (int) ( ( alpha_rad - aa_min[ind_sizex] ) / ( ( aa_max[ind_sizex] - aa_min[ind_sizex] ) / 10.0 ) );
0452
0453 h_xres_npix_beta_alpha[ind_sizex][ind_beta][ind_alpha]->Fill( rechitresx );
0454
0455
0456
0457
0458
0459
0460
0461
0462
0463 }
0464
0465 }
0466 else if ( subdetId == 2 )
0467 {
0468
0469 int sizey = nypix;
0470
0471 if ( !( sizey == 1 && bigy == 1 ) )
0472 {
0473 if ( sizey > 2 ) sizey = 2;
0474
0475 int ind_sizey = sizey - 1;
0476 int ind_beta = -9999;
0477
0478 if ( betap_rad < 0.3 ) ind_beta = 0;
0479 else if ( betap_rad > 0.4 ) ind_beta = 9;
0480 else
0481 ind_beta = (int) ( ( betap_rad - 0.3 ) / ( ( 0.4 - 0.3 ) / 10.0 ) );
0482
0483 h_forward_yres_npix_beta[ind_sizey][ind_beta]->Fill( rechitresy );
0484
0485 }
0486
0487
0488 int sizex = nxpix;
0489
0490 if ( !( sizex == 1 && bigx == 1 ) )
0491 {
0492 if ( sizex > 2 ) sizex = 2;
0493
0494 int ind_sizex = sizex - 1;
0495 int ind_alpha = -9999;
0496
0497 if ( sizex == 1 )
0498 {
0499 if ( alphap_rad < 0.15 ) ind_alpha = 0;
0500 else if ( alphap_rad > 0.30 ) ind_alpha = 9;
0501 else
0502 ind_alpha = (int) ( ( alphap_rad - 0.15 ) / ( ( 0.3 - 0.15 ) / 10.0 ) );
0503 }
0504 if ( sizex > 1 )
0505 {
0506 if ( alphap_rad < 0.15 ) ind_alpha = 0;
0507 else if ( alphap_rad > 0.50 ) ind_alpha = 9;
0508 else
0509 ind_alpha = (int) ( ( alphap_rad - 0.15 ) / ( ( 0.5 - 0.15 ) / 10.0 ) );
0510 }
0511
0512 h_forward_xres_npix_alpha[ind_sizex][ind_alpha]->Fill( rechitresx );
0513
0514 }
0515
0516 }
0517 else
0518 cout << " Wrong Detector ID !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! " << endl;
0519
0520 }
0521
0522 float low = -999.9;
0523 float high = -999.9;
0524 if ( do_residuals )
0525 {
0526 low = -0.02;
0527 high = 0.02;
0528 }
0529 else
0530 {
0531 low = -10.0;
0532 high = 10.0;
0533 }
0534 TF1* myfunc = new TF1("myfunc", ng, low, high, 3);
0535 myfunc->SetParNames("norm", "mean", "sigma");
0536 myfunc->SetParameter(0, 100.0);
0537 myfunc->FixParameter(1, 0.0);
0538 myfunc->SetParameter(2, 0.00020);
0539 myfunc->SetLineColor(kRed);
0540
0541 Double_t sigma = -99999.9;
0542 Double_t ssigma = -99999.9;
0543
0544
0545
0546 if ( do_yb )
0547 {
0548 for (int k=0; k<6; ++k)
0549 {
0550 for (int i=0; i<4; ++i)
0551 for (int j=0; j<10; ++j)
0552 {
0553 can_y_barrel_sizey_alpha[k][i]->cd(j + 1);
0554 h_yres_npix_alpha_beta[k][i][j]->Draw();
0555
0556 double n = h_yres_npix_alpha_beta[k][i][j]->GetEntries();
0557 double rms = h_yres_npix_alpha_beta[k][i][j]->GetRMS();
0558 double binw = h_yres_npix_alpha_beta[k][i][j]->GetBinWidth(1);
0559
0560 if ( n != 0.0 )
0561 {
0562 double norm = n*binw;
0563
0564 myfunc->SetParameter(0, norm);
0565 myfunc->FixParameter(1, 0.0);
0566 myfunc->SetParameter(2, rms/3.0);
0567
0568
0569
0570 if ( k==0 && j<3 )
0571 h_yres_npix_alpha_beta[k][i][j]->Fit("myfunc", "LQR");
0572 else
0573 h_yres_npix_alpha_beta[k][i][j]->Fit("myfunc", "QR");
0574
0575 sigma = myfunc->GetParameter(2);
0576 ssigma = myfunc->GetParError(2);
0577 h_yres_npix_alpha[k][i]->SetBinContent(j+1, sigma);
0578 h_yres_npix_alpha[k][i]->SetBinError(j+1, ssigma);
0579
0580 h_yres_npix_alpha_rms[k][i]->SetBinContent(j+1, rms);
0581 }
0582
0583 if ( sigma < 0.0 )
0584 {
0585 sigma = rms;
0586 cout << "Bad error, check fit convergence !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!" << endl;
0587 }
0588
0589 fprintf( datfile,
0590 "%d %d %d %d %f %f \n",
0591 1, k, i, j, sigma, rms );
0592
0593 }
0594
0595 for (int i=0; i<4; ++i)
0596 {
0597 if ( !do_residuals )
0598 {
0599 h_yres_npix_alpha[k][i]->SetMinimum(0.0);
0600 h_yres_npix_alpha[k][i]->SetMinimum(2.0);
0601 }
0602 can_y_barrel_sizey[k]->cd(i+1);
0603 h_yres_npix_alpha_rms[k][i]->Draw("");
0604 h_yres_npix_alpha[k][i]->Draw("same");
0605 }
0606 }
0607
0608 if ( do_plots )
0609 {
0610 for (int i=0; i<6; ++i)
0611 {
0612 for (int j=0; j<4; ++j)
0613 {
0614 if ( do_residuals )
0615 sprintf(hname, "res_y_barrel_sizey_alpha_%i_%i.eps", i, j);
0616 else
0617 sprintf(hname, "pull_y_barrel_sizey_alpha_%i_%i.eps", i, j);
0618
0619 can_y_barrel_sizey_alpha[i][j]->SaveAs(hname);
0620 }
0621
0622 if ( do_residuals )
0623 sprintf(hname, "res_y_barrel_sizey_%i.eps", i);
0624 else
0625 sprintf(hname, "pull_y_barrel_sizey_%i.eps", i);
0626
0627 can_y_barrel_sizey[i]->SaveAs(hname);
0628 }
0629 }
0630 }
0631
0632
0633
0634
0635
0636 if ( do_xb )
0637 {
0638 for (int k=0; k<3; ++k)
0639 {
0640 for (int i=0; i<4; ++i)
0641 for (int j=0; j<10; ++j)
0642 {
0643 Double_t sigma = -99999.9;
0644 can_x_barrel_sizex_beta[k][i]->cd(j + 1);
0645 h_xres_npix_beta_alpha[k][i][j]->Draw();
0646
0647 double n = h_xres_npix_beta_alpha[k][i][j]->GetEntries();
0648 double rms = h_xres_npix_beta_alpha[k][i][j]->GetRMS();
0649 double binw = h_xres_npix_beta_alpha[k][i][j]->GetBinWidth(1);
0650
0651 if ( n != 0.0 )
0652 {
0653 double norm = n*binw;
0654
0655 myfunc->SetParameter(0, norm);
0656 myfunc->FixParameter(1, 0.0);
0657 myfunc->SetParameter(2, rms/3.0);
0658
0659 h_xres_npix_beta_alpha[k][i][j]->Fit("myfunc", "QR");
0660 sigma = myfunc->GetParameter(2);
0661 ssigma = myfunc->GetParError(2);
0662 h_xres_npix_beta[k][i]->SetBinContent(j+1, sigma);
0663 h_xres_npix_beta[k][i]->SetBinError(j+1, ssigma);
0664
0665 h_xres_npix_beta_rms[k][i]->SetBinContent(j+1, rms);
0666 }
0667
0668 if ( sigma < 0.0 )
0669 {
0670 sigma = rms;
0671 cout << "Bad error, check fit convergence !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!" << endl;
0672 }
0673
0674 fprintf( datfile,
0675 "%d %d %d %d %f %f \n",
0676 2, k, i, j, sigma, rms );
0677
0678 }
0679
0680 for (int i=0; i<4; ++i)
0681 {
0682 if ( !do_residuals )
0683 {
0684 h_xres_npix_beta[k][i]->SetMinimum(0.0);
0685 h_xres_npix_beta[k][i]->SetMinimum(2.0);
0686 }
0687 can_x_barrel_sizex[k]->cd(i+1);
0688 h_xres_npix_beta_rms[k][i]->Draw("");
0689 h_xres_npix_beta[k][i]->Draw("same");
0690 }
0691 }
0692
0693 if ( do_plots )
0694 {
0695 for (int i=0; i<3; ++i)
0696 {
0697 for (int j=0; j<4; ++j)
0698 {
0699 if ( do_residuals )
0700 sprintf(hname, "res_x_barrel_sizex_beta_%i_%i.eps", i, j);
0701 else
0702 sprintf(hname, "pull_x_barrel_sizex_beta_%i_%i.eps", i, j);
0703
0704 can_x_barrel_sizex_beta[i][j]->SaveAs(hname);
0705 }
0706
0707 if ( do_residuals )
0708 sprintf(hname, "res_x_barrel_sizex_%i.eps", i);
0709 else
0710 sprintf(hname, "pull_x_barrel_sizex_%i.eps", i);
0711
0712 can_x_barrel_sizex[i]->SaveAs(hname);
0713 }
0714 }
0715
0716 }
0717
0718
0719
0720
0721 if ( do_yf )
0722 {
0723 for ( int k=0; k<2; ++k )
0724 for (int i=0; i<10; ++i)
0725 {
0726 can_y_forward_sizey[k]->cd(i+1);
0727 h_forward_yres_npix_beta[k][i]->Draw();
0728
0729 double n = h_forward_yres_npix_beta[k][i]->GetEntries();
0730 double rms = h_forward_yres_npix_beta[k][i]->GetRMS();
0731 double binw = h_forward_yres_npix_beta[k][i]->GetBinWidth(1);
0732
0733 if ( n != 0.0 )
0734 {
0735 double norm = n*binw;
0736
0737 myfunc->SetParameter(0, norm);
0738 myfunc->FixParameter(1, 0.0);
0739 myfunc->SetParameter(2, rms/3.0);
0740
0741 h_forward_yres_npix_beta[k][i]->Fit("myfunc", "QR");
0742 sigma = myfunc->GetParameter(2);
0743 ssigma = myfunc->GetParError(2);
0744 h_forward_yres_npix[k]->SetBinContent(i+1, sigma);
0745 h_forward_yres_npix[k]->SetBinError(i+1, ssigma);
0746
0747 h_forward_yres_npix_rms[k]->SetBinContent(i+1, rms);
0748 }
0749
0750 if ( sigma < 0.0 )
0751 {
0752 sigma = rms;
0753 cout << "Bad error, check fit convergence !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!" << endl;
0754 }
0755
0756 fprintf( datfile,
0757 "%d %d %d %d %f %f \n",
0758 3, 0, k, i, sigma, rms );
0759
0760 }
0761
0762 for (int i=0; i<2; ++i)
0763 {
0764 if ( !do_residuals )
0765 {
0766 h_forward_yres_npix[i]->SetMinimum(0.0);
0767 h_forward_yres_npix[i]->SetMinimum(2.0);
0768 }
0769 can_y_forward->cd(i+1);
0770 h_forward_yres_npix_rms[i]->Draw("");
0771 h_forward_yres_npix[i]->Draw("same");
0772 }
0773
0774 if ( do_plots )
0775 {
0776 for (int i=0; i<2; ++i)
0777 {
0778 if ( do_residuals )
0779 sprintf(hname, "res_y_forward_sizey_%i.eps", i);
0780 else
0781 sprintf(hname, "pull_y_forward_sizey_%i.eps", i);
0782
0783 can_y_forward_sizey[i]->SaveAs(hname);
0784 }
0785
0786 if ( do_residuals )
0787 sprintf(hname, "res_y_forward.eps");
0788 else
0789 sprintf(hname, "pull_y_forward.eps");
0790
0791 can_y_forward->SaveAs(hname);
0792 }
0793 }
0794
0795
0796
0797
0798 if ( do_xf )
0799 {
0800 for (int k=0; k<2; ++k)
0801 for (int i=0; i<10; ++i)
0802 {
0803 can_x_forward_sizex[k]->cd(i+1);
0804 h_forward_xres_npix_alpha[k][i]->Draw();
0805
0806 double n = h_forward_xres_npix_alpha[k][i]->GetEntries();
0807 double rms = h_forward_xres_npix_alpha[k][i]->GetRMS();
0808 double binw = h_forward_xres_npix_alpha[k][i]->GetBinWidth(1);
0809
0810 if ( n != 0.0 )
0811 {
0812 double norm = n*binw;
0813
0814 myfunc->SetParameter(0, norm);
0815 myfunc->FixParameter(1, 0.0);
0816 myfunc->SetParameter(2, rms/3.0);
0817
0818 h_forward_xres_npix_alpha[k][i]->Fit("myfunc", "QR");
0819 sigma = myfunc->GetParameter(2);
0820 ssigma = myfunc->GetParError(2);
0821 h_forward_xres_npix[k]->SetBinContent(i+1, sigma);
0822 h_forward_xres_npix[k]->SetBinError(i+1, ssigma);
0823
0824 h_forward_xres_npix_rms[k]->SetBinContent(i+1, rms);
0825
0826 }
0827
0828 if ( sigma < 0.0 )
0829 {
0830 sigma = rms;
0831 cout << "Bad error, check fit convergence !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!" << endl;
0832 }
0833
0834 fprintf( datfile,
0835 "%d %d %d %d %f %f \n",
0836 4, 0, k, i, sigma, rms );
0837
0838 }
0839
0840 for (int i=0; i<2; ++i)
0841 {
0842 if ( !do_residuals )
0843 {
0844 h_forward_xres_npix[i]->SetMinimum(0.0);
0845 h_forward_xres_npix[i]->SetMinimum(2.0);
0846 }
0847 can_x_forward->cd(i+1);
0848 h_forward_xres_npix_rms[i]->Draw("");
0849 h_forward_xres_npix[i]->Draw("same");
0850 }
0851
0852 if ( do_plots )
0853 {
0854 for (int i=0; i<2; ++i)
0855 {
0856 if ( do_residuals )
0857 sprintf(hname, "res_x_forward_sizex_%i.eps", i);
0858 else
0859 sprintf(hname, "pull_x_forward_sizex_%i.eps", i);
0860
0861 can_x_forward_sizex[i]->SaveAs(hname);
0862 }
0863
0864 if ( do_residuals )
0865 sprintf(hname, "res_x_forward.eps");
0866 else
0867 sprintf(hname, "pull_x_forward.eps");
0868
0869 can_x_forward->SaveAs(hname);
0870 }
0871
0872 }
0873
0874 }