EMCal/EMCALmacros/InterlopationV4.C

0001 #include <cmath>
0002 #include <TFile.h>
0003 #include <TString.h>
0004 #include <TLine.h>
0005 #include <TTree.h>
0006 #include <TLatex.h>
0007 #include <TGraphErrors.h>
0008 #include <cassert>
0009 #include <iostream>
0010 #include <fstream>
0011 using namespace std;
0012
0013
0014
0015 using std::cout;
0016 using std::endl;
0017 #endif
0018
0019 void InterlopationV4()
0020 {
0021
0022      gStyle->SetOptFit(0);
0023    gStyle->SetOptStat(0);
0024
0025     TFile *fin = new TFile("EnergyPosition2.root");
0026
0027 //Set Tower Horizontal Lines//
0028 double y1 = 90;
0029 double y2 = 113;
0030 double y3 = 141;
0031 double y4 = 166;
0032 double y5 = 190;
0033
0034 //Set Tower Vertical Lines//
0035
0036 double x1 = 204;
0037 double x2 = 227.5;
0038 double x3 = 248;
0039 double x4 = 265;
0040 double x5 = 283;
0041
0042
0043
0044 int binx;
0045
0046 int forward = 2;
0047 int backward = 2;
0048
0049
0050 int biny;
0051 double step = 1;
0052
0053 double Xmin = 170;
0054 double Xmax = 310.0;
0055
0056 double Ymin = 60;
0057 double Ymax = 245.0;
0058
0059
0060
0061
0062 cout << Ymax << endl;
0063
0064 int NumX = (Xmax - Xmin)/step;
0065 int NumY = (Ymax - Ymin)/step;
0066
0067 const int XBins = NumX;
0068
0069 const int YBins = NumY;
0070 double x;
0071 double y;
0072
0073 cout << "XBins = " << XBins << endl;
0074 cout << "YBins = " << YBins << endl;
0075
0076
0077 double value[XBins][YBins];
0078
0079 double finalvalue[XBins][YBins];
0080
0081 TH2D *Inter = new TH2D("Inter","",XBins,Xmin,Xmax,YBins,Ymin,Ymax);
0082
0083 /*
0084 double ylower = 35;
0085 double yhigher = 188;
0086 double xlower = 178;
0087 double xhigher = 301;
0088 */
0089 int k;
0090
0091 int NumberofValue;
0092 double Sum;
0093
0094
0095
0096 TCanvas *c1 = new TCanvas("c1", "c1",0,0,800,600);
0097
0098 c1->cd();
0099
0100 EnPo->GetXaxis()->SetTitle("Horizontal Axis (mm)");
0101 EnPo->GetYaxis()->SetTitle("Vertical Axis (mm)");
0102
0103 EnPo->SetTitle("Energy vs Horizontal and Vertical Positions Before Interpolation - 10 Degree");
0104
0105 EnPo->Draw("colz");
0106 c1->Update();
0107
0108 c1->SaveAs("Before10Degree.png");
0109
0110
0111 for(int i = 0; i < XBins; i++)
0112 {
0113
0114 //x = Xmin + step*i;
0115
0116
0117 for(int j = 0; j < YBins; j++)
0118 {
0119
0120 //  TCanvas *c1 = new TCanvas("c1", "c1",0,0,800,600);
0121 //  y = Ymin + step*j;
0122
0123
0124
0125
0126 value[i][j] = EnPo->GetBinContent(i,j);
0127
0128
0129 }
0130 }
0131
0132 //Horizontal Interpolation Begins//
0133
0134 for(int i = backward; i < XBins - backward; i++)
0135 {
0136
0137     for(int j = backward; j < YBins - backward; j++)
0138 {
0139
0140 //cout << "OK" << endl;
0141
0142 NumberofValue = 0;
0143 Sum = 0;
0144 for(int k = i - backward; k < i + forward +1; k++)
0145 {
0146
0147     for(int m = j - backward; m < j+backward+1; m++)
0148     {
0149
0150         if(value[k][m] > 0) NumberofValue = NumberofValue + 1;
0151
0152         Sum = Sum + value[k][m];
0153
0154             }
0155
0156 }
0157
0158 cout << "Number of Value = " << NumberofValue << endl;
0159
0160 cout << "Sum = " << Sum << endl;
0161
0162 if(NumberofValue > 0) finalvalue[i][j] = Sum/NumberofValue;
0163
0164 cout << " i = " << i <<  "  j = " << j << endl;
0165
0166 cout << "Final Value is " << finalvalue[i][j]  << endl;
0167
0168
0169 Inter->SetBinContent(i,j,finalvalue[i][j]);
0170
0171
0172 }
0173
0174
0175 }
0176
0177 /*
0178 cout << "Horizontal Done" << endl;
0179
0180
0181 //Vertical Interpolation Begins//
0182
0183 cout << "DONE Value" << endl;
0184
0185 for(int i = 1; i < XBins; i++)
0186 {
0187
0188     for(int j = 1; j < YBins; j++)
0189 {
0190
0191     cout << "Before Updating: Binx =  " << i << "   " << "Bin y = " << j << " Value =  " << value[i][j] << endl;
0192
0193
0194 if (value[i][j] < 6.8 && value[i][j-1] >0 && i < XBins && j < YBins)
0195 {
0196
0197
0198 k = j;
0199 do {
0200     k = k + 1;
0201
0202 //  cout << "k value is = " << k << endl;
0203 }
0204 while(value[i][k] == 0 && k < YBins);
0205
0206 //cout << "neW k " << k  << endl;
0207
0208 value[i][j] =value[i][j-1] + (value[i][k] - value[i][j-1])/(k+1-j);
0209
0210     cout << "After Updating: Binx =  " << i << "   " << "  Bin y = " << j << "  Value =  " <<value[i][j] << endl;
0211
0212 }
0213
0214
0215 }
0216
0217
0218 }
0219
0220 cout << "Vk =" << k << endl;
0221
0222 //Vertical Interpolation Ends//
0223
0224 cout << "VERTICAL DONE" << endl;
0225
0226 */
0227 /*
0228 for(int i =0; i < XBins; i++)
0229 {
0230
0231
0232
0233     for(int j = 1; j < YBins; j++)
0234 {
0235
0236
0237 }
0238
0239
0240 }
0241 */
0242 Inter->GetXaxis()->SetTitle("Horizontal Axis (mm)");
0243 Inter->GetYaxis()->SetTitle("Vertical Axis (mm)");
0244 Inter->SetTitle("Energy vs Horizontal and Vertical Positions Before Interpolation - 10 Degree");
0245
0246   TFile *fout = new TFile("Interpolated2.root","RECREATE");
0247 Inter->Write();
0248
0249
0250 TCanvas *c22 = new TCanvas("c22", "c22",0,0,800,600);
0251
0252 c22->cd();
0253
0254
0255
0256 Inter->SetMinimum(4.5);
0257
0258 Inter->SetMaximum(9);
0259
0260 Inter->Draw("colz");
0261
0262 // Horizontal line Tower 1 //
0263    TLine *l1 = new TLine(Xmin,y1,Xmax,y1);
0264    l1->SetLineWidth(1);
0265    l1->Draw("same");
0266
0267
0268  TLatex* texCol;
0269 texCol= new TLatex(Xmin,y1,"Tower 6");
0270   texCol->SetTextAlign(13);
0271   texCol->SetTextSize(0.03);
0272  texCol->SetTextFont(42);
0273   texCol->Draw();
0274
0275
0276 // Horizontal line Tower 2 //
0277
0278
0279    TLine *l2 = new TLine(Xmin,y2,Xmax,y2);
0280    l2->SetLineWidth(1);
0281    l2->Draw("same");
0282
0283
0284
0285 texCol= new TLatex(Xmin,y2,"Tower 5");
0286   texCol->SetTextAlign(13);
0287   texCol->SetTextSize(0.03);
0288  texCol->SetTextFont(42);
0289   texCol->Draw();
0290
0291 // Horizontal line Tower 3//
0292
0293
0294    TLine *l3 = new TLine(Xmin,y3,Xmax,y3);
0295    l3->SetLineWidth(1);
0296    l3->Draw("same");
0297
0298
0299
0300 texCol= new TLatex(Xmin,y3,"Tower 4");
0301   texCol->SetTextAlign(13);
0302   texCol->SetTextSize(0.03);
0303  texCol->SetTextFont(42);
0304   texCol->Draw();
0305
0306
0307
0308 // Horizontal line Tower 4 //
0309
0310
0311    TLine *l4 = new TLine(Xmin,y4,Xmax,y4);
0312    l4->SetLineWidth(1);
0313    l4->Draw("same");
0314
0315
0316
0317 texCol= new TLatex(Xmin,y4,"Tower 3");
0318   texCol->SetTextAlign(13);
0319   texCol->SetTextSize(0.03);
0320  texCol->SetTextFont(42);
0321   texCol->Draw();
0322
0323   // Horizontal line Tower 5 //
0324
0325
0326    TLine *l4 = new TLine(Xmin,y5,Xmax,y5);
0327    l4->SetLineWidth(1);
0328    l4->Draw("same");
0329
0330
0331
0332
0333 texCol= new TLatex(Xmin,y5,"Tower 2");
0334   texCol->SetTextAlign(13);
0335   texCol->SetTextSize(0.03);
0336  texCol->SetTextFont(42);
0337   texCol->Draw();
0338
0339
0340   // Horizontal line Tower 6 //
0341
0342      TLine *l12 = new TLine(Xmax,Ymax,Xmax,Ymax);
0343    l12->SetLineWidth(1);
0344    l12->Draw("same");
0345
0346
0347  TLatex* texCol;
0348 texCol= new TLatex(Xmin,215,"Tower 1");
0349   texCol->SetTextAlign(13);
0350   texCol->SetTextSize(0.03);
0351  texCol->SetTextFont(42);
0352   texCol->Draw();
0353
0354
0355
0356
0357
0358 // Vertical line Tower 1 //
0359    TLine *l6 = new TLine(x1,Ymin,x1,Ymax);
0360    l6->SetLineWidth(1);
0361    l6->Draw("same");
0362
0363
0364  TLatex* texCol;
0365 texCol= new TLatex(x1,Ymax,"Tower 6");
0366   texCol->SetTextAlign(33);
0367   texCol->SetTextSize(0.03);
0368  texCol->SetTextFont(42);
0369   texCol->Draw();
0370
0371
0372 // Horizontal line Tower 2 //
0373
0374
0375    TLine *l7 = new TLine(x2,Ymin,x2,Ymax);
0376    l7->SetLineWidth(1);
0377    l7->Draw("same");
0378
0379
0380  TLatex* texCol;
0381 texCol= new TLatex(x2,Ymax,"Tower 5");
0382   texCol->SetTextAlign(33);
0383   texCol->SetTextSize(0.03);
0384  texCol->SetTextFont(42);
0385   texCol->Draw();
0386
0387 // Horizontal line Tower 3//
0388
0389
0390    TLine *l8 = new TLine(x3,Ymin,x3,Ymax);
0391    l8->SetLineWidth(1);
0392    l8->Draw("same");
0393
0394
0395  TLatex* texCol;
0396 texCol= new TLatex(x3,Ymax,"Tower 4");
0397   texCol->SetTextAlign(33);
0398   texCol->SetTextSize(0.03);
0399  texCol->SetTextFont(42);
0400   texCol->Draw();
0401
0402
0403 // Horizontal line Tower 4 //
0404
0405    TLine *l9 = new TLine(x4,Ymin,x4,Ymax);
0406    l9->SetLineWidth(1);
0407    l9->Draw("same");
0408
0409
0410  TLatex* texCol;
0411 texCol= new TLatex(x4,Ymax,"Tower 3");
0412   texCol->SetTextAlign(33);
0413   texCol->SetTextSize(0.03);
0414  texCol->SetTextFont(42);
0415   texCol->Draw();
0416
0417
0418   // Horizontal line Tower 5 //
0419
0420
0421    TLine *l10 = new TLine(x5,Ymin,x5,Ymax);
0422    l10->SetLineWidth(1);
0423    l10->Draw("same");
0424
0425
0426  TLatex* texCol;
0427 texCol= new TLatex(x5,Ymax,"Tower 2");
0428   texCol->SetTextAlign(33);
0429   texCol->SetTextSize(0.03);
0430  texCol->SetTextFont(42);
0431   texCol->Draw();
0432
0433   // Horizontal line Tower 6 //
0434
0435      TLine *l11 = new TLine(Xmax,Ymin,Xmax,Ymax);
0436    l11->SetLineWidth(1);
0437    l11->Draw("same");
0438
0439
0440  TLatex* texCol;
0441 texCol= new TLatex(Xmax,Ymax,"Tower 1");
0442   texCol->SetTextAlign(33);
0443   texCol->SetTextSize(0.03);
0444  texCol->SetTextFont(42);
0445   texCol->Draw();
0446
0447
0448
0449 c22->Update();
0450
0451 c22->SaveAs("Interpolation10Degree.png");
0452
0453 }