Changeset 9953 for trunk/Cosy/tpoint/plot_m1.C
- Timestamp:
- 09/30/10 09:57:16 (14 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/Cosy/tpoint/plot_m1.C
r9512 r9953 1 #include <iomanip>1 #include "mars/MTime.h" 2 2 3 void DrawMarker(TVirtualPad *pad, Double_t r0, Double_t phi0, Double_t r1, Double_t phi1) 4 { 5 TView *view = pad->GetView(); 3 using namespace std; 6 4 7 if (!view)8 {9 cout << "No View!" << endl;10 return;11 }12 13 TMarker mark0;14 TMarker mark1;15 mark0.SetMarkerStyle(kStar);16 mark1.SetMarkerStyle(kStar);17 mark1.SetMarkerColor(kRed);18 19 r0 /= 90;20 r1 /= 90;21 phi0 *= TMath::DegToRad();22 phi1 *= TMath::DegToRad();23 24 Double_t x0[3] = { r0*cos(phi0), r0*sin(phi0), 0};25 Double_t x1[3] = { r1*cos(phi1), r1*sin(phi1), 0};26 27 mark0.DrawMarker(x0[0], x0[1]);28 mark1.DrawMarker(x1[0], x1[1]);29 30 return;31 Double_t y0[3], y1[3];32 33 view->WCtoNDC(x0, y0);34 view->WCtoNDC(x1, y1);35 36 mark0.DrawMarker(y0[0], y0[1]);37 mark1.DrawMarker(y1[0], y1[1]);38 }39 40 int fill(const char *fname, TGraph *g, TH1 *h)41 {42 5 /* 43 TH2F h2res1("Res2D1", " Dataset positions on the sky ", 32, 0, 360, 10, 0, 90);44 h2res1.SetBit(TH1::kNoStats);45 h2res1.DrawCopy("surf1pol");46 gPad->Modified();47 gPad->Update();48 gPad->SetTheta(90);49 gPad->SetPhi(-90);50 51 DrawMarker(gPad, 45, 0, 0, 0);52 gPad->Modified();53 gPad->Update();54 55 return;56 57 */58 ifstream fin(fname);59 60 cout << "Reading " << setw(23) << gSystem->BaseName(fname) << "..." << flush;61 62 while (1)63 {64 TString str;65 str.ReadLine(fin);66 if (!fin)67 break;68 69 if (str.Contains("#"))70 continue;71 72 Float_t alt, az, dalt, daz, mjd;73 sscanf(str.Data(), "%f %f %*f %*f %*f %*f %f %f %f",74 &az, &alt, &dalt, &daz, &mjd);75 76 if (dalt==0/* || GetResidual(alt, az, alt+dalt, az+daz)>0.1*/)77 continue;78 79 mjd -= 53140.097505;80 81 Double_t res = MAstro::GetDevAbs(90-(alt-dalt), -dalt, -daz);82 83 g[0].SetPoint(g[0].GetN(), g[0].GetN(), fabs(dalt));84 g[1].SetPoint(g[1].GetN(), g[1].GetN(), fabs(daz));85 g[2].SetPoint(g[2].GetN(), g[2].GetN(), res);86 87 h->Fill(res);88 }89 90 cout << "done (" << setw(3) << (Int_t)h->GetEntries() << "/";91 cout << setw(4) << g[0].GetN() << ") " << flush;92 93 return g[0].GetN();94 }95 96 struct Description_t97 {98 const char *fName;99 const char *fTitle;100 const char *fFile;101 };102 103 const Int_t counts = 29+10+18+1+13+11+31-27+3+5+18;104 Description_t desc[counts] =105 {106 /*107 // 29. Apr. 2004 ~25800108 // 5. Aug. 2004 ~32000109 // 19. Aug. 2004 ~33530110 // 7. Jun. 2005 ~57650111 // 8. Jun. 2005112 // 9. Jun. 2005 ~57860113 // 12. Sep. 2005 ~68338114 // 24. Nov. 2005 ~75562115 // 17. Oct. 2006 ~103130116 // 17. Jun. 2007 ~248193117 */118 119 6 // Culmination tests 120 7 //{"0411", "TPoints Residuals 11/2004" , "tpoint/tpoint0411.txt"}, … … 146 33 //{"05111", "TPoints Residuals 11/2005-1" , "tpoint/tpoint0511-1.txt"}, 147 34 148 // New pointing model installed (24.11.2005) 149 {"05112", "TPoints Residuals 11/2005-2" , "tpoint/tpoint0511-2.txt"}, 150 {"+0512", "TPoints Residuals 12/2005" , "tpoint/tpoint0512.txt"}, 151 {"+0601", "TPoints Residuals 1/2006" , "tpoint/tpoint0601.txt"}, 152 {"+0603-1", "TPoints Residuals 3/2006-1" , "tpoint/tpoint0603-1.txt"}, 35 // [...] 36 */ 153 37 154 // Changes to the mirror 155 {"0603-2", "TPoints Residuals 3/2006-2" , "tpoint/tpoint0603-2.txt"}, 156 {"+0604", "TPoints Residuals 4/2006" , "tpoint/tpoint0604.txt"}, 157 {"+0607", "TPoints Residuals 7/2006" , "tpoint/tpoint0607.txt"}, 158 {"+0608", "TPoints Residuals 8/2006" , "tpoint/tpoint0608.txt"}, 159 {"+0609", "TPoints Residuals 9/2006" , "tpoint/tpoint0609.txt"}, 160 {"+0610", "TPoints Residuals 10/2006" , "tpoint/tpoint0610.txt"}, 38 Double_t dates[] = { 39 MTime(2005, 3, 20).GetAxisTime(), 40 MTime(2005, 4, 29).GetAxisTime(), 41 MTime(2005, 5, 25).GetAxisTime(), 42 MTime(2005, 6, 8).GetAxisTime(), // New pointing model 43 MTime(2005, 8, 15).GetAxisTime(), 44 // MTime(2005, 9, 12).GetAxisTime(), // New pointing model 45 MTime(2005, 11, 10).GetAxisTime(), // New mirror alignment after Tenerife meeting 46 MTime(2005, 11, 24).GetAxisTime(), // New pointing model 47 MTime(2006, 3, 19).GetAxisTime(), // Changes to the mirror 48 // 2006, 4, 23 // Mirror refocussing 49 MTime(2006, 10, 17).GetAxisTime(), // New pointing model 50 MTime(2007, 6, 17).GetAxisTime(), // New pointing model 51 MTime(2007, 8, 4).GetAxisTime(), // Mirror refocussing 52 MTime(2007, 10, 18).GetAxisTime(), // New pointing model 53 MTime(2008, 1, 14).GetAxisTime(), // New pointing model 54 MTime(2008, 6, 11).GetAxisTime(), // New pointing model 55 MTime(2008, 6, 18).GetAxisTime(), // New pointing model 56 // 2009, 3, 7 // New pointing model 57 // 2009, 5, 14 // New pointing model 58 // Are we missing TPoints between 080618 and 090501?? 59 MTime(2009, 5, 1).GetAxisTime(), // Drive upgrade started 60 MTime(2009, 5, 11).GetAxisTime(), // Upgrade finished 61 MTime(2009, 5, 12).GetAxisTime(), // First new pointing model 62 MTime(2009, 5, 13).GetAxisTime(), // Second new pointing model 63 MTime(2009, 6, 11).GetAxisTime(), // Jump (reason unknown) 64 MTime(2009, 7, 23).GetAxisTime(), // New LUTs 65 MTime(2009, 8, 17).GetAxisTime(), // New pointing model 66 MTime(2010, 02, 01).GetAxisTime(), // New pointing model 67 // 2010, 02, 03 // New LUTs for M2(!) 68 MTime(2010, 02, 27).GetAxisTime(), // New pointing model 69 MTime(2010, 03, 31).GetAxisTime(), // New pointing model 70 // 2010, 06, 14 // New LUTs 71 MTime(2010, 07, 03).GetAxisTime(), // New LUTs 72 MTime(2010, 9, 29).GetAxisTime(), // New pointing model 73 MTime(-1).GetAxisTime(), 74 -1 75 }; 161 76 162 // New pointing model: 6/10/17 163 {"0611", "TPoints Residuals 11/2006" , "tpoint/tpoint0611.txt"}, 164 {"+0612", "TPoints Residuals 12/2006" , "tpoint/tpoint0612.txt"}, 165 {"+0701", "TPoints Residuals 1/2007" , "tpoint/tpoint0701.txt"}, 166 {"+0702", "TPoints Residuals 2/2007" , "tpoint/tpoint0702.txt"}, 167 {"+0703", "TPoints Residuals 3/2007" , "tpoint/tpoint0703.txt"}, 168 169 {"+0704", "TPoints Residuals 4/2007" , "tpoint/tpoint0704.txt"}, 170 {"+0705", "TPoints Residuals 5/2007" , "tpoint/tpoint0705.txt"}, 171 {"+0706", "TPoints Residuals 6/2007" , "tpoint/tpoint0706.txt"}, 172 {"+07071", "TPoints Residuals 7/2007-1" , "tpoint/tpoint0707-1.txt"}, 173 {"+07072", "TPoints Residuals 7/2007-2" , "tpoint/tpoint0707-2.txt"}, 174 175 // New pointing model: 07/06/17 176 {"", "", ""}, 177 178 // AMC adjust: 07/08/04 179 {"0708", "TPoints Residuals 7/8/23", "tpoint/m1/tpoints20070823.txt"}, 180 181 // New pointing model: 07/10/18 182 {"", "", ""}, 183 184 // New pointing model: 08/01/14 185 {"0801", "TPoints Residuals 7/8/23", "tpoint/m1/tpoints20080115.txt"}, 186 187 // New pointing model 080611 188 {"0806", "TPoints Residuals 7/8/23", "tpoint/tpoint0806.txt"}, 189 190 // New pointing model 080618 191 {"0807", "", "tpoint/tpoint0807.txt"}, 192 {"+0807", "", "tpoint/tpoint0808.txt"}, 193 {"+0807", "", "tpoint/tpoint0811.txt"}, 194 {"+0807", "", "tpoint/tpoint0812.txt"}, 195 {"+0901", "", "tpoint/tpoint0901.txt"}, 196 197 // Drive upgrade started 198 {"", "", ""}, 199 200 // Drive Upgrade finished 090511 201 {"090512", "TPoints Residuals 7/8/23", "tpoint/m1/2009_05_12/tpoint_20090511_220118.txt"}, 202 {"+090512", "TPoints Residuals 7/8/23", "tpoint/m1/2009_05_12/tpoint_20090511_221333.txt"}, 203 {"+090512", "TPoints Residuals 7/8/23", "tpoint/m1/2009_05_12/tpoint_20090511_222434.txt"}, 204 {"+090512", "TPoints Residuals 7/8/23", "tpoint/m1/2009_05_12/tpoint_20090511_223152.txt"}, 205 {"+090512", "TPoints Residuals 7/8/23", "tpoint/m1/2009_05_12/tpoint_20090511_223706.txt"}, 206 {"+090512", "TPoints Residuals 7/8/23", "tpoint/m1/2009_05_12/tpoint_20090511_224056.txt"}, 207 {"+090512", "TPoints Residuals 7/8/23", "tpoint/m1/2009_05_12/tpoint_20090511_231933.txt"}, 208 {"+090512", "TPoints Residuals 7/8/23", "tpoint/m1/2009_05_12/tpoint_20090512_011131.txt"}, 209 210 // First new pointing model 090512 211 {"090513", "TPoints Residuals 7/8/23", "tpoint/m1/2009_05_13/tpoint_20090512_210644.txt"}, 212 {"+090513", "TPoints Residuals 7/8/23", "tpoint/m1/2009_05_13/tpoint_20090513_002757.txt"}, 213 {"+090513", "TPoints Residuals 7/8/23", "tpoint/m1/2009_05_13/tpoint_20090513_025124.txt"}, 214 {"+090513", "TPoints Residuals 7/8/23", "tpoint/m1/2009_05_13/tpoint_20090513_033505.txt"}, 215 216 // Second new pointing model 090513 217 {"090514", "TPoints 09/05/14", "tpoint/m1/2009_05_14/tpoint_20090514_011903.txt"}, 218 219 {"090517", "TPoints 09/05/17" ,"tpoint/m1/2009_05_17/tpoint_20090516_234825.txt"}, 220 {"+090517","TPoints 09/05/17", "tpoint/m1/2009_05_17/tpoint_20090517_023340.txt"}, 221 // New period 222 {"+090517","TPoints 09/05/17", "tpoint/m1/2009_06_01/tpoint_20090531_215148.txt"}, 223 {"+090517","TPoints 09/05/17", "tpoint/m1/2009_06_01/tpoint_20090531_222549.txt"}, 224 {"+090517","TPoints 09/05/17", "tpoint/m1/2009_06_02/tpoint_20090601_223009.txt"}, 225 {"+090517","TPoints 09/05/17", "tpoint/m1/2009_06_03/tpoint_20090602_213509.txt"}, 226 {"+090517","TPoints 09/05/17", "tpoint/m1/2009_06_03/tpoint_20090603_011936.txt"}, 227 {"+090517","TPoints 09/05/17", "tpoint/m1/2009_06_04/tpoint_20090603_215840.txt"}, 228 {"+090517","TPoints 09/05/17", "tpoint/m1/2009_06_04/tpoint_20090603_230510.txt"}, 229 {"+090517","TPoints 09/05/17", "tpoint/m1/2009_06_05/tpoint_20090604_215943.txt"}, 230 {"+090517","TPoints 09/05/17", "tpoint/m1/2009_06_05/tpoint_20090604_232320.txt"}, 231 {"+090517","TPoints 09/05/17", "tpoint/m1/2009_06_10/tpoint_20090609_231948.txt"}, 232 {"+090517","TPoints 09/05/17", "tpoint/m1/2009_06_11/tpoint_20090611_011155.txt"}, 233 234 // ------- Something happened (reason unknown) -------- 235 236 {"090613","TPoints 09/05/17", "tpoint/m1/2009_06_11/tpoint_20090611_023625.txt"}, 237 238 {"+090613","TPoints 09/05/17", "tpoint/m1/2009_06_13/tpoint_20090613_033838.txt"}, 239 {"+090613","TPoints 09/05/17", "tpoint/m1/2009_06_14/tpoint_20090614_021257.txt"}, 240 {"+090613","TPoints 09/05/17", "tpoint/m1/2009_06_18/tpoint_20090618_041433.txt"}, 241 242 {"+090613","TPoints 09/05/17", "tpoint/m1/2009_07_02/tpoint_20090701_215304.txt"}, 243 {"+090613","TPoints 09/05/17", "tpoint/m1/2009_07_02/tpoint_20090701_222059.txt"}, 244 {"+090613","TPoints 09/05/17", "tpoint/m1/2009_07_02/tpoint_20090701_224051.txt"}, 245 {"+090613","TPoints 09/05/17", "tpoint/m1/2009_07_02/tpoint_20090701_225615.txt"}, 246 {"+090613","TPoints 09/05/17", "tpoint/m1/2009_07_02/tpoint_20090701_230946.txt"}, 247 248 {"+090613","TPoints 09/05/17", "tpoint/m1/2009_07_03/tpoint_20090702_225940.txt"}, 249 {"+090613","TPoints 09/05/17", "tpoint/m1/2009_07_04/tpoint_20090703_224721.txt"}, 250 {"+090613","TPoints 09/05/17", "tpoint/m1/2009_07_05/tpoint_20090705_012638.txt"}, 251 252 {"+090613","TPoints 09/05/17", "tpoint/m1/2009_07_12/tpoint_20090712_014300.txt"}, 253 {"+090613","TPoints 09/05/17", "tpoint/m1/2009_07_12/tpoint_20090712_024710.txt"}, 254 {"+090613","TPoints 09/05/17", "tpoint/m1/2009_07_13/tpoint_20090713_004241.txt"}, 255 {"+090613","TPoints 09/05/17", "tpoint/m1/2009_07_13/tpoint_20090713_025934.txt"}, 256 257 {"+090613","TPoints 09/05/17", "tpoint/m1/2009_07_14/tpoint_20090714_024729.txt"}, 258 {"+090613","TPoints 09/05/17", "tpoint/m1/2009_07_15/tpoint_20090715_021320.txt"}, 259 {"+090613","TPoints 09/05/17", "tpoint/m1/2009_07_15/tpoint_20090715_025237.txt"}, 260 261 // ----- AMC adjust (23.7.) ----- 262 263 {"090723","TPoints 09/05/17", "tpoint/m1/2009_07_24/tpoint_20090724_051142.txt"}, 264 {"+090723", "TPoints 08/2009" , "tpoint/m1/2009_08_02/tpoint_20090801_223848.txt"}, 265 {"+090723","TPoints 08/2009" , "tpoint/m1/2009_08_03/tpoint_20090802_224434.txt"}, 266 {"+090723","TPoints 08/2009" , "tpoint/m1/2009_08_04/tpoint_20090803_230214.txt"}, 267 {"+090723","TPoints 08/2009" , "tpoint/m1/2009_08_10/tpoint_20090809_233524.txt"}, 268 {"+090723","TPoints 08/2009" , "tpoint/m1/2009_08_10/tpoint_20090810_014642.txt"}, 269 {"+090723","TPoints 08/2009" , "tpoint/m1/2009_08_10/tpoint_20090810_020616.txt"}, 270 {"+090723","TPoints 08/2009" , "tpoint/m1/2009_08_11/tpoint_20090810_232543.txt"}, 271 {"+090723","TPoints 08/2009" , "tpoint/m1/2009_08_12/tpoint_20090812_001846.txt"}, 272 {"+090723","TPoints 08/2009" , "tpoint/m1/2009_08_13/tpoint_20090813_030807.txt"}, 273 274 // ------ 09/08/17 new pointing model ------- 275 276 {"090817","TPoints", "tpoint/m1/2009_09_01/tpoint_20090831_233034.txt"}, 277 {"+090817","TPoints", "tpoint/m1/2009_09_02/tpoint_20090901_231918.txt"}, 278 {"+090817","TPoints", "tpoint/m1/2009_09_02/tpoint_20090902_005204.txt"}, 279 {"+090817","TPoints", "tpoint/m1/2009_09_08/tpoint_20090907_213404.txt"}, 280 {"+090817","TPoints", "tpoint/m1/2009_09_09/tpoint_20090908_225027.txt"}, 281 {"+090817","TPoints", "tpoint/m1/2009_09_09/tpoint_20090908_231720.txt"}, 282 {"+090817","TPoints", "tpoint/m1/2009_09_09/tpoint_20090909_010649.txt"}, 283 {"+090817","TPoints", "tpoint/m1/2009_09_10/tpoint_20090909_235222.txt"}, 284 {"+090817","TPoints", "tpoint/m1/2009_09_10/tpoint_20090910_025054.txt"}, 285 {"+090817","TPoints", "tpoint/m1/2009_09_12/tpoint_20090912_004826.txt"}, 286 {"+090817","TPoints", "tpoint/m1/2009_09_17/tpoint_20090917_015404.txt"}, 287 {"+090817","TPoints", "tpoint/m1/2009_09_18/tpoint_20090917_202449.txt"}, 288 {"+090817","TPoints", "tpoint/m1/2009_09_23/tpoint_20090922_223445.txt"}, 289 {"+090817","TPoints", "tpoint/m1/2009_09_24/tpoint_20090923_204444.txt"}, 290 {"+090817","TPoints", "tpoint/m1/2009_09_24/tpoint_20090924_041824.txt"}, 291 {"+090817","TPoints", "tpoint/m1/2009_09_25/tpoint_20090924_225725.txt"}, 292 {"+090817","TPoints", "tpoint/m1/2009_09_26/tpoint_20090925_212537.txt"}, 293 {"+090817","TPoints", "tpoint/m1/2009_09_27/tpoint_20090926_212332.txt"}, 294 }; 77 #include "plot.C" 295 78 296 79 void plot_m1() 297 80 { 298 TGraph g[3];299 81 300 MBinning bins(100, 0, 0.2); 82 MDirIter Next; 83 Next.AddDirectory("tpoint/m1", "tpoint*.txt", -1); 301 84 302 TArrayI n(counts); 303 304 TH1F hx[counts]; 305 Int_t num = -1; 306 for (int i=0; i<counts; i++) 307 { 308 if (desc[i].fName[0]!='+') 309 { 310 num++; 311 312 hx[num].SetNameTitle(desc[i].fName, desc[i].fTitle); 313 hx[num].SetDirectory(0); 314 bins.Apply(hx[num]); 315 } 316 317 cout << setw(2) << num << ": " << flush; 318 n[num] = fill(desc[i].fFile, g, &hx[num]); 319 cout << " Mean: " << setw(5) << setprecision(3) << hx[num].GetMean()*3600 << " sec +/- " << hx[num].GetRMS()*3600 << endl; 320 } 321 322 n.Set(++num); 323 324 g[0].SetMarkerColor(kGreen); 325 g[1].SetMarkerColor(kMagenta); 326 g[2].SetMarkerColor(kBlack); 327 //g[2].SetLineColor(kBlack); 328 g[0].SetMarkerStyle(kFullDotMedium); 329 g[1].SetMarkerStyle(kFullDotMedium); 330 g[2].SetMarkerStyle(kFullDotLarge); 331 332 // --------- First Canvas ---------- 333 334 new TCanvas("Vs Time", ""); 335 336 TObject *obj[4]; 337 338 for (int i=0; i<3; i++) 339 { 340 g[i].SetFillColor(kWhite); 341 g[i].SetLineColor(kWhite); 342 obj[i] = g[i].Clone(); 343 obj[i]->SetBit(kCanDelete); 344 obj[i]->Draw(i==0?"AP":"P"); 345 } 346 347 TLegend leg(0.905, 0.86, 0.99, 0.99); 348 leg.AddEntry(obj[0], " \\Delta\\theta"); 349 leg.AddEntry(obj[1], " \\Delta\\phi"); 350 leg.AddEntry(obj[2], " \\Delta"); 351 leg.DrawClone()->SetBit(kCanDelete); 352 353 TLine l; 354 l.SetLineColor(kBlue); 355 for (int i=0; i<n.GetSize(); i++) 356 l.DrawLine(n[i], 0, n[i], 0.2); 357 358 // --------- Second Canvas ---------- 359 360 new TCanvas("Distrib", ""); 361 362 Double_t max=0; 363 for (int i=0; i<n.GetSize(); i++) 364 { 365 if (hx[i].GetEntries()==0) 366 { 367 cout << "Skip #" << i << endl; 368 continue; 369 } 370 371 hx[i].Scale(1./hx[i].GetEntries()); 372 max = TMath::Max(max, hx[i].GetMaximum()); 373 } 374 for (int i=0; i<n.GetSize(); i++) 375 { 376 hx[i].SetMaximum(max*1.05); 377 if (i<6) 378 hx[i].SetLineColor(kRed+i); 379 else 380 hx[i].SetMarkerStyle(kPlus+i-6); 381 382 383 } 384 385 for (int i=0; i<counts; i++) 386 hx[i].DrawCopy(i==0?"LP":"LPsame"); 387 388 Double_t time[] = { 389 MTime(2005, 3, 20).GetAxisTime(), 390 MTime(2005, 4, 29).GetAxisTime(), 391 MTime(2005, 5, 25).GetAxisTime(), 392 MTime(2005, 6, 8).GetAxisTime(), 393 MTime(2005, 8, 15).GetAxisTime(), 394 // MTime(2005, 9, 12).GetAxisTime(), 395 MTime(2005, 11, 10).GetAxisTime(), 396 MTime(2005, 11, 24).GetAxisTime(), 397 MTime(2006, 3, 19).GetAxisTime(), 398 MTime(2006, 10, 17).GetAxisTime(), 399 MTime(2007, 6, 17).GetAxisTime(), 400 MTime(2007, 8, 4).GetAxisTime(), 401 MTime(2007, 10, 18).GetAxisTime(), 402 MTime(2008, 1, 14).GetAxisTime(), 403 MTime(2008, 6, 11).GetAxisTime(), 404 MTime(2008, 6, 18).GetAxisTime(), 405 // Are we missing TPoints between 080618 and 090501?? 406 MTime(2009, 5, 1).GetAxisTime(), 407 MTime(2009, 5, 11).GetAxisTime(), 408 MTime(2009, 5, 12).GetAxisTime(), 409 MTime(2009, 5, 13).GetAxisTime(), 410 MTime(2009, 5, 16).GetAxisTime(), 411 MTime(2009, 6, 11).GetAxisTime(), 412 MTime(2009, 7, 23).GetAxisTime(), 413 MTime(2009, 8, 17).GetAxisTime(), 414 MTime(2009, 12, 31).GetAxisTime(), 415 }; 416 417 TH1D histres[4]; 418 419 MBinning bins; 420 bins.SetEdges(TArrayD(24, time)); 421 422 bins.Apply(histres[0]); 423 bins.Apply(histres[1]); 424 bins.Apply(histres[2]); 425 bins.Apply(histres[3]); 426 427 TGraphAsymmErrors result[4]; 428 TGraph resultm; 429 for (int i=0; i<n.GetSize(); i++) 430 { 431 cout << i+1 << " - Mean: " << Form("%5.1f +- %5.1f", hx[i].GetMean()*3600, hx[i].GetRMS()*3600); 432 cout << " (Overflows=" << hx[i].GetBinContent(hx[i].GetNbinsX()+1)*hx[i].GetEntries() << ") " << (int)hx[i].GetEntries() << endl; 433 434 if (hx[i].GetEntries()<1) 435 continue; 436 437 /* 438 TF1 fg("fg", "gaus", 0, 1); 439 hx[i].Fit(&fg); 440 result2.SetPoint(result.GetN(), result.GetN()+7, hx[i].GetMean()); 441 result2.SetPointError(result.GetN()-1, 0, hx[i].GetRMS()/sqrt(hx[i].GetEntries())); 442 result.SetPoint(result.GetN(), result.GetN()+7, fg.GetParameter(1)); 443 result.SetPointError(result.GetN()-1, 0, fg.GetParameter(2)); 444 */ 445 446 //Double_t q[4] = { MMath::GaussProb(0.5), MMath::GaussProb(1), MMath::GaussProb(2), MMath::GaussProb(3) }; 447 Double_t q[4] = { 0.5, MMath::GaussProb(1), MMath::GaussProb(2), MMath::GaussProb(3) }; 448 Double_t rc[4]; 449 hx[i].GetQuantiles(4, rc, q); 450 451 for (int j=0; j<4; j++) 452 { 453 histres[j].SetBinContent(i+1, 60*rc[j]); 454 455 result[j].SetPoint(i, time[i], 60*rc[j]); 456 result[j].SetPointError(i, 0, 0, 60*rc[j], 0); 457 458 //result[j].SetPoint(result[j].GetN(), result[j].GetN()+1, 60*rc[j]); 459 //result[j].SetPointError(result[j].GetN()-1, 0, 0, 60*rc[j], 0); 460 } 461 462 // result.SetPoint(result.GetN(), result.GetN()+1, rc[1]); 463 // result2.SetPoint(result2.GetN(), result2.GetN()+1, rc[2]); 464 465 //result.SetPointError(result.GetN()-1, 0.5, 0);//rc[2]-rc[1]); 466 //result2.SetPointError(result.GetN()-1, 0.5, 0);//rc[2]-rc[1]); 467 468 // result.SetPointError(result.GetN()-1, 0, 0, 0/*rc[1]-rc[0]*/, rc[2]-rc[1]); 469 // result2.SetPointError(result.GetN()-1, 0, 0, 0/*rc[1]-rc[0]*/, rc[2]-rc[1]); 470 471 // result2.SetPointError(result.GetN()-1, 0, 0, rc[2], 0); 472 473 // resultm.SetPoint(resultm.GetN(), resultm.GetN()+1, 60*hx[i].GetMean()); 474 resultm.SetPoint(resultm.GetN(), (time[i]+time[i+1])/2, 60*hx[i].GetMean()); 475 476 } 477 478 new TCanvas; 479 480 gPad->SetBorderMode(0); 481 gPad->SetFrameBorderMode(0); 482 gPad->SetFillColor(kWhite); 483 gPad->SetRightMargin(0.01); 484 gPad->SetTopMargin(0.01); 485 gPad->SetLeftMargin(0.06); 486 gPad->SetGridy(); 487 488 //Int_t col[] = { 12, 15, 17, 19 }; 489 //Int_t col[] = { 12, 16, 18, 0 }; 490 Int_t col[] = { 13, 16, 19, 0 }; 491 492 TH1 *h = &histres[2];//result[3].GetHistogram(); 493 h->SetXTitle(""); 494 h->SetYTitle("Residual / arcmin"); 495 h->SetBit(TH1::kNoStats); 496 h->GetXaxis()->CenterTitle(); 497 h->GetYaxis()->CenterTitle(); 498 h->GetYaxis()->SetTitleOffset(0.75); 499 // h->GetXaxis()->SetTimeFormat("%m/%y %F1995-01-01 00:00:00 GMT"); 500 // h->GetXaxis()->SetTimeDisplay(1); 501 502 h->GetXaxis()->SetLabelColor(kWhite); 503 504 TLine line; 505 506 for (int j=2; j>=0; j--) 507 { 508 histres[j].SetMinimum(0); 509 histres[j].SetFillColor(col[j]);//12+2*j); 510 511 histres[j].DrawCopy(j==2?"":"same"); 512 513 /* 514 //result[j].SetLineColor(kBlue); 515 //result[j].SetLineWidth(2); 516 //result[j].SetMarkerColor(kBlue); 517 result[j].SetMinimum(0); 518 //result2.SetMarkerStyle(kFullDotMedium); 519 //result[j].SetMarkerStyle(23); 520 result[j].SetFillColor(col[j]);//12+2*j); 521 result[j].DrawClone(j==3 ? "ABX" : "B"); // E3 B 522 */ 523 } 524 525 526 resultm.SetMarkerStyle(20); 527 resultm.SetMarkerSize(0.8); 528 resultm.DrawClone("P"); 529 530 line.DrawLine(time[0], 0, time[0], histres[2].GetMaximum()); 531 for (int i=0; i<bins.GetNumBins()-1; i++) 532 line.DrawLine(time[i+1], 0, time[i+1], histres[2].GetBinContent(i+1)); 533 534 TText txt; 535 txt.SetTextSize(0.037); 536 // txt.SetTextAngle(-45); 537 txt.SetTextAlign(23); 538 for (int m=4; m<13; m++) 539 { 540 Double_t monl = MTime(2005, m, 1,0).GetAxisTime(); 541 Double_t mont = MTime(2005, m, 15,0).GetAxisTime(); 542 // txt.DrawText(mon, -0.12, Form("%02d/05", m)); 543 txt.DrawText(mont, -0.10, Form("%d", m)); 544 line.DrawLine(monl, -0.12, monl, 0.12); 545 } 546 for (int m=1; m<13; m++) 547 { 548 Double_t monl = MTime(2006, m, 1,0).GetAxisTime(); 549 Double_t mont = MTime(2006, m, 15,0).GetAxisTime(); 550 // txt.DrawText(mon, -0.12, Form("%02d/06", m)); 551 txt.DrawText(mont, -0.10, Form("%d", m)); 552 line.DrawLine(monl, -0.12, monl, 0.12); 553 } 554 for (int m=1; m<13; m++) 555 { 556 Double_t monl = MTime(2007, m, 1,0).GetAxisTime(); 557 Double_t mont = MTime(2007, m, 15,0).GetAxisTime(); 558 // txt.DrawText(mon, -0.12, Form("%02d/07", m)); 559 txt.DrawText(mont, -0.10, Form("%d", m)); 560 line.DrawLine(monl, -0.12, monl, 0.12); 561 } 562 for (int m=1; m<13; m++) 563 { 564 Double_t monl = MTime(2008, m, 1,0).GetAxisTime(); 565 Double_t mont = MTime(2008, m, 15,0).GetAxisTime(); 566 // txt.DrawText(mon, -0.12, Form("%02d/07", m)); 567 txt.DrawText(mont, -0.10, Form("%d", m)); 568 line.DrawLine(monl, -0.12, monl, 0.12); 569 } 570 571 for (int m=1; m<13; m++) 572 { 573 Double_t monl = MTime(2009, m, 1,0).GetAxisTime(); 574 Double_t mont = MTime(2009, m, 15,0).GetAxisTime(); 575 // txt.DrawText(mon, -0.12, Form("%02d/07", m)); 576 txt.DrawText(mont, -0.10, Form("%d", m)); 577 line.DrawLine(monl, -0.12, monl, 0.12); 578 } 579 580 Double_t y6 = MTime(2006,1,1,0).GetAxisTime(); 581 Double_t y7 = MTime(2007,1,1,0).GetAxisTime(); 582 Double_t y8 = MTime(2008,1,1,0).GetAxisTime(); 583 Double_t y9 = MTime(2009,1,1,0).GetAxisTime(); 584 Double_t y0 = MTime(2010,1,1,0).GetAxisTime(); 585 // Double_t y0 = MTime(2009,5,1,0).GetAxisTime(); 586 587 txt.SetTextSize(0.042); 588 txt.DrawText(y6-(y7-y6)/2, -0.6, "2005"); 589 txt.DrawText((y6+y7)/2, -0.6, "2006"); 590 txt.DrawText((y7+y8)/2, -0.6, "2007"); 591 txt.DrawText((y8+y9)/2, -0.6, "2008"); 592 txt.DrawText((y9+y0)/2, -0.6, "2009"); 593 594 line.DrawLine(y6, -0.7, y6, 0.26); 595 line.DrawLine(y7, -0.7, y7, 0.26); 596 line.DrawLine(y8, -0.7, y8, 0.26); 597 line.DrawLine(y9, -0.7, y9, 0.26); 598 line.SetLineStyle(3); 599 line.DrawLine(y6, 0.26, y6, 1.05*histres[2].GetMaximum()); 600 line.DrawLine(y7, 0.26, y7, 1.05*histres[2].GetMaximum()); 601 line.DrawLine(y8, 0.26, y8, 1.05*histres[2].GetMaximum()); 602 line.DrawLine(y9, 0.26, y9, 1.05*histres[2].GetMaximum()); 603 604 line.SetLineColor(kBlue); 605 line.SetLineWidth(2); 606 line.SetLineStyle(kSolid); 607 line.DrawLine(time[0], 1*360/16384.*60, time[bins.GetNumBins()], 1*360/16384.*60); 608 line.SetLineStyle(9); 609 line.DrawLine(time[0], 2*360/16384.*60, time[bins.GetNumBins()], 2*360/16384.*60); 610 line.SetLineStyle(7); 611 line.DrawLine(time[0], 3*360/16384.*60, time[bins.GetNumBins()], 3*360/16384.*60); 612 613 /* result.SetMinimum(-0.06); 614 result.SetMarkerStyle(kFullDotMedium); 615 result.DrawClone("LP");*/ 85 plot(Next); 616 86 }
Note:
See TracChangeset
for help on using the changeset viewer.