Ignore:
Timestamp:
09/30/10 09:57:16 (14 years ago)
Author:
tbretz
Message:
Simplified and unified plotting macros.
File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/Cosy/tpoint/plot_m1.C

    r9512 r9953  
    1 #include <iomanip>
     1#include "mars/MTime.h"
    22
    3 void DrawMarker(TVirtualPad *pad, Double_t r0, Double_t phi0, Double_t r1, Double_t phi1)
    4     {
    5         TView *view = pad->GetView();
     3using namespace std;
    64
    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 {
    425/*
    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_t
    97 {
    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    ~25800
    108      //    5. Aug. 2004    ~32000
    109      //   19. Aug. 2004    ~33530
    110      //    7. Jun. 2005    ~57650
    111      //    8. Jun. 2005
    112      //    9. Jun. 2005    ~57860
    113      //   12. Sep. 2005    ~68338
    114      //   24. Nov. 2005    ~75562
    115      //   17. Oct. 2006   ~103130
    116      //   17. Jun. 2007   ~248193
    117      */
    118 
    1196    // Culmination tests
    1207    //{"0411",  "TPoints Residuals 11/2004" ,  "tpoint/tpoint0411.txt"},
     
    14633    //{"05111", "TPoints Residuals 11/2005-1" , "tpoint/tpoint0511-1.txt"},
    14734
    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*/
    15337
    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"},
     38Double_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};
    16176
    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"
    29578
    29679void plot_m1()
    29780{
    298     TGraph g[3];
    29981
    300     MBinning bins(100, 0, 0.2);
     82    MDirIter Next;
     83    Next.AddDirectory("tpoint/m1", "tpoint*.txt", -1);
    30184
    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);
    61686}
Note: See TracChangeset for help on using the changeset viewer.