Changeset 3929 for trunk/MagicSoft/Mars/mjobs
- Timestamp:
- 05/03/04 10:11:15 (21 years ago)
- Location:
- trunk/MagicSoft/Mars/mjobs
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/MagicSoft/Mars/mjobs/MGCamDisplays.cc
r3889 r3929 39 39 #include "MGeomCam.h" 40 40 #include "TVirtualPad.h" 41 #include "TProfile.h"42 #include "TF1.h"43 41 44 42 #include "MLog.h" … … 99 97 gPad->SetBorderMode(0); 100 98 gPad->SetTicks(); 101 DrawRadialProfile(obj1);99 obj1->DrawRadialProfile(); 102 100 } 103 101 … … 109 107 gPad->SetBorderMode(0); 110 108 gPad->SetTicks(); 111 DrawProjection(obj1,fit);109 obj1->DrawProjection(fit); 112 110 } 113 111 114 // --------------------------------------------------------------------------115 //116 // Draw a projection of MHCamera onto the y-axis values. Depending on the117 // variable fit, the following fits are performed:118 //119 // 1: Single Gauss (for distributions flat-fielded over the whole camera)120 // 2: Double Gauss (for distributions different for inner and outer pixels)121 // 3: Triple Gauss (for distributions with inner, outer pixels and outliers)122 // 4: flat (for the probability distributions)123 // (1-4:) Moreover, sectors 6,1 and 2 of the camera and sectors 3,4 and 5 are124 // drawn separately, for inner and outer pixels.125 // 5: Fit Inner and Outer pixels separately by a single Gaussian126 // (only for MAGIC cameras)127 // 6: Fit Inner and Outer pixels separately by a single Gaussian and display128 // additionally the two camera halfs separately (for MAGIC camera)129 //130 //131 void MGCamDisplays::DrawProjection(MHCamera *obj, Int_t fit) const132 {133 134 TArrayI inner(1);135 inner[0] = 0;136 137 TArrayI outer(1);138 outer[0] = 1;139 140 if (fit==5 || fit==6)141 {142 143 if (obj->GetGeomCam().InheritsFrom("MGeomCamMagic"))144 {145 TArrayI s0(6);146 s0[0] = 6;147 s0[1] = 1;148 s0[2] = 2;149 s0[3] = 3;150 s0[4] = 4;151 s0[5] = 5;152 112 153 TArrayI s1(3);154 s1[0] = 6;155 s1[1] = 1;156 s1[2] = 2;157 158 TArrayI s2(3);159 s2[0] = 3;160 s2[1] = 4;161 s2[2] = 5;162 163 gPad->Clear();164 TVirtualPad *pad = gPad;165 pad->Divide(2,1);166 167 TH1D *inout[2];168 inout[0] = obj->ProjectionS(s0, inner, "Inner");169 inout[1] = obj->ProjectionS(s0, outer, "Outer");170 171 inout[0]->SetDirectory(NULL);172 inout[1]->SetDirectory(NULL);173 174 for (int i=0; i<2; i++)175 {176 pad->cd(i+1);177 inout[i]->SetLineColor(kRed+i);178 inout[i]->SetBit(kCanDelete);179 inout[i]->Draw();180 inout[i]->Fit("gaus","Q");181 182 if (fit == 6)183 {184 TH1D *half[2];185 half[0] = obj->ProjectionS(s1, i==0 ? inner : outer , "Sector 6-1-2");186 half[1] = obj->ProjectionS(s2, i==0 ? inner : outer , "Sector 3-4-5");187 188 for (int j=0; j<2; j++)189 {190 half[j]->SetLineColor(kRed+i+j);191 half[j]->SetDirectory(0);192 half[j]->SetBit(kCanDelete);193 half[j]->Draw("same");194 }195 }196 197 }198 199 gLog << all << obj->GetName()200 << Form("%s%5.3f%s%3.2f"," Mean: Inner Pixels: ",201 inout[0]->GetFunction("gaus")->GetParameter(1),"+-",202 inout[0]->GetFunction("gaus")->GetParError(1));203 gLog << Form("%s%5.3f%s%3.2f"," Outer Pixels: ",204 inout[1]->GetFunction("gaus")->GetParameter(1),"+-",205 inout[1]->GetFunction("gaus")->GetParError(1)) << endl;206 gLog << all << obj->GetName()207 << Form("%s%5.3f%s%3.2f"," Sigma: Inner Pixels: ",208 inout[0]->GetFunction("gaus")->GetParameter(2),"+-",209 inout[0]->GetFunction("gaus")->GetParError(2));210 gLog << Form("%s%5.3f%s%3.2f"," Outer Pixels: ",211 inout[1]->GetFunction("gaus")->GetParameter(2),"+-",212 inout[1]->GetFunction("gaus")->GetParError(2)) << endl;213 214 }215 return;216 }217 218 TH1D *obj2 = (TH1D*)obj->Projection(obj->GetName());219 obj2->SetDirectory(0);220 obj2->Draw();221 obj2->SetBit(kCanDelete);222 223 224 if (obj->GetGeomCam().InheritsFrom("MGeomCamMagic"))225 {226 TArrayI s0(3);227 s0[0] = 6;228 s0[1] = 1;229 s0[2] = 2;230 231 TArrayI s1(3);232 s1[0] = 3;233 s1[1] = 4;234 s1[2] = 5;235 236 237 TH1D *halfInOut[4];238 239 // Just to get the right (maximum) binning240 halfInOut[0] = obj->ProjectionS(s0, inner, "Sector 6-1-2 Inner");241 halfInOut[1] = obj->ProjectionS(s1, inner, "Sector 3-4-5 Inner");242 halfInOut[2] = obj->ProjectionS(s0, outer, "Sector 6-1-2 Outer");243 halfInOut[3] = obj->ProjectionS(s1, outer, "Sector 3-4-5 Outer");244 245 for (int i=0; i<4; i++)246 {247 halfInOut[i]->SetLineColor(kRed+i);248 halfInOut[i]->SetDirectory(0);249 halfInOut[i]->SetBit(kCanDelete);250 halfInOut[i]->Draw("same");251 }252 }253 254 const Double_t min = obj2->GetBinCenter(obj2->GetXaxis()->GetFirst());255 const Double_t max = obj2->GetBinCenter(obj2->GetXaxis()->GetLast());256 const Double_t integ = obj2->Integral("width")/2.5;257 const Double_t mean = obj2->GetMean();258 const Double_t rms = obj2->GetRMS();259 const Double_t width = max-min;260 261 const TString dgausformula = "([0]-[3])/[2]*exp(-0.5*(x-[1])*(x-[1])/[2]/[2])"262 "+[3]/[5]*exp(-0.5*(x-[4])*(x-[4])/[5]/[5])";263 264 const TString tgausformula = "([0]-[3]-[6])/[2]*exp(-0.5*(x-[1])*(x-[1])/[2]/[2])"265 "+[3]/[5]*exp(-0.5*(x-[4])*(x-[4])/[5]/[5])"266 "+[6]/[8]*exp(-0.5*(x-[7])*(x-[7])/[8]/[8])";267 TF1 *f=0;268 switch (fit)269 {270 case 1:271 f = new TF1("sgaus", "gaus(0)", min, max);272 f->SetLineColor(kYellow);273 f->SetBit(kCanDelete);274 f->SetParNames("Area", "#mu", "#sigma");275 f->SetParameters(integ/rms, mean, rms);276 f->SetParLimits(0, 0, integ);277 f->SetParLimits(1, min, max);278 f->SetParLimits(2, 0, width/1.5);279 280 obj2->Fit(f, "QLR");281 break;282 283 case 2:284 f = new TF1("dgaus",dgausformula.Data(),min,max);285 f->SetLineColor(kYellow);286 f->SetBit(kCanDelete);287 f->SetParNames("A_{tot}", "#mu1", "#sigma1", "A2", "#mu2", "#sigma2");288 f->SetParameters(integ,(min+mean)/2.,width/4.,289 integ/width/2.,(max+mean)/2.,width/4.);290 // The left-sided Gauss291 f->SetParLimits(0,integ-1.5 , integ+1.5);292 f->SetParLimits(1,min+(width/10.), mean);293 f->SetParLimits(2,0 , width/2.);294 // The right-sided Gauss295 f->SetParLimits(3,0 , integ);296 f->SetParLimits(4,mean, max-(width/10.));297 f->SetParLimits(5,0 , width/2.);298 obj2->Fit(f,"QLRM");299 break;300 301 case 3:302 f = new TF1("tgaus",tgausformula.Data(),min,max);303 f->SetLineColor(kYellow);304 f->SetBit(kCanDelete);305 f->SetParNames("A_{tot}","#mu_{1}","#sigma_{1}",306 "A_{2}","#mu_{2}","#sigma_{2}",307 "A_{3}","#mu_{3}","#sigma_{3}");308 f->SetParameters(integ,(min+mean)/2,width/4.,309 integ/width/3.,(max+mean)/2.,width/4.,310 integ/width/3.,mean,width/2.);311 // The left-sided Gauss312 f->SetParLimits(0,integ-1.5,integ+1.5);313 f->SetParLimits(1,min+(width/10.),mean);314 f->SetParLimits(2,width/15.,width/2.);315 // The right-sided Gauss316 f->SetParLimits(3,0.,integ);317 f->SetParLimits(4,mean,max-(width/10.));318 f->SetParLimits(5,width/15.,width/2.);319 // The Gauss describing the outliers320 f->SetParLimits(6,0.,integ);321 f->SetParLimits(7,min,max);322 f->SetParLimits(8,width/4.,width/1.5);323 obj2->Fit(f,"QLRM");324 break;325 326 case 4:327 obj2->Fit("pol0", "Q");328 obj2->GetFunction("pol0")->SetLineColor(kYellow);329 break;330 331 case 9:332 break;333 334 default:335 obj2->Fit("gaus", "Q");336 obj2->GetFunction("gaus")->SetLineColor(kYellow);337 break;338 }339 }340 341 // --------------------------------------------------------------------------342 //343 // Draw a projection of MHCamera vs. the radius from the central pixel.344 //345 // The inner and outer pixels are drawn separately, both fitted by a polynomial346 // of grade 1.347 //348 void MGCamDisplays::DrawRadialProfile(MHCamera *obj) const349 {350 351 TProfile *obj2 = (TProfile*)obj->RadialProfile(obj->GetName());352 obj2->SetDirectory(0);353 obj2->Draw();354 obj2->SetBit(kCanDelete);355 356 if (obj->GetGeomCam().InheritsFrom("MGeomCamMagic"))357 {358 359 TArrayI s0(6);360 s0[0] = 1;361 s0[1] = 2;362 s0[2] = 3;363 s0[3] = 4;364 s0[4] = 5;365 s0[5] = 6;366 367 TArrayI inner(1);368 inner[0] = 0;369 370 TArrayI outer(1);371 outer[0] = 1;372 373 // Just to get the right (maximum) binning374 TProfile *half[2];375 half[0] = obj->RadialProfileS(s0, inner,Form("%s%s",obj->GetName(),"Inner"));376 half[1] = obj->RadialProfileS(s0, outer,Form("%s%s",obj->GetName(),"Outer"));377 378 for (Int_t i=0; i<2; i++)379 {380 Double_t min = obj->GetGeomCam().GetMinRadius(i);381 Double_t max = obj->GetGeomCam().GetMaxRadius(i);382 383 half[i]->SetLineColor(kRed+i);384 half[i]->SetDirectory(0);385 half[i]->SetBit(kCanDelete);386 half[i]->Draw("same");387 half[i]->Fit("pol1","Q","",min,max);388 half[i]->GetFunction("pol1")->SetLineColor(kRed+i);389 half[i]->GetFunction("pol1")->SetLineWidth(1);390 }391 }392 }393 -
trunk/MagicSoft/Mars/mjobs/MGCamDisplays.h
r3875 r3929 8 8 class TCanvas; 9 9 class MHCamera; 10 class MCamEvent;11 10 class MGCamDisplays 12 11 { 13 12 protected: 14 13 15 void DrawProjection ( MHCamera *obj, Int_t fit) const; // Draw projection of pixels values16 void DrawRadialProfile( MHCamera *obj) const; // Draw projection of pixels values onto camera radius17 14 void CamDraw(TCanvas &c, const Int_t x, const Int_t y, const MHCamera &cam1, 18 15 const Int_t fit, const Int_t rad=0,
Note:
See TracChangeset
for help on using the changeset viewer.