Ignore:
Timestamp:
05/04/04 18:22:12 (21 years ago)
Author:
rico
Message:
*** empty log message ***
File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/MagicSoft/Mars/mtemp/mifae/programs/makeHillas.cc

    r3947 r3965  
    1919#include "MHCalibrationChargeCam.h"
    2020#include "MHCalibrationRelTimeCam.h"
     21#include "MExtractor.h"
     22#include "MExtractFixedWindow.h"
    2123#include "MExtractSignal.h"
    2224#include "MCalibrationChargeCalc.h"
    … …  
    5355
    5456// initial and final time slices to be used in signal extraction
    55 const Byte_t hifirst = 1;
    56 const Byte_t hilast  = 14;
    57 const Byte_t lofirst = 4;
    58 const Byte_t lolast  = 13;
     57const Byte_t hifirst = 0;
     58const Byte_t hilast  = 13;
     59const Byte_t lofirst = 3;
     60const Byte_t lolast  = 12;
    5961
    6062// declaration of variables read from datacards
    … …  
    7173
    7274const TString defaultcard="input.datacard";
     75/*************************************************************/
    7376static void Usage()
    7477{
    … …  
    7982  gLog << "     -?/-h: This help" << endl << endl;
    8083}
    81 
     84/*************************************************************/
    8285int main(int argc, char **argv)
    8386{
    … …  
    105108void makeHillas()
    106109{
    107 #if 0
    108   MStatusDisplay *display = new MStatusDisplay;
    109   display->SetUpdateTime(3000);
    110   display->Resize(850,700);
    111  
    112   gStyle->SetOptStat(1111);
    113   gStyle->SetOptFit();
    114 #endif
    115  
     110  // Set the general tasks/containers
     111  MExtractFixedWindow    extractor;
     112  extractor.SetRange(hifirst,hilast,lofirst,lolast);
     113
     114  MCalibrationQECam   qecam;
     115  MBadPixelsCam       badcam;
     116  MGeomCamMagic       geomcam;
     117  MGeomApply          geomapl;
     118
    116119  /************************************/
    117120  /* FIRST LOOP: PEDESTAL COMPUTATION */
    118121  /************************************/
    119  
    120   MParList plist1;
    121   MTaskList tlist1;
    122   plist1.AddToList(&tlist1);
    123  
    124   // containers
    125   MPedestalCam   pedcam;
    126   MBadPixelsCam  badcam; 
    127   // badcam.AsciiRead("badpixels.dat");
    128 
    129   plist1.AddToList(&pedcam);
    130   plist1.AddToList(&badcam);
    131 
    132   //tasks
    133   MReadMarsFile read("Events");
    134   static_cast<MRead&>(read).AddFiles(pediter);
    135   read.DisableAutoScheme();
    136  
    137   MGeomApply     geomapl;
    138   MPedCalcPedRun pedcalc;
    139 
    140   tlist1.AddToList(&read);
    141   tlist1.AddToList(&geomapl);
    142   tlist1.AddToList(&pedcalc);
    143 
    144   // Create and setup the eventloop
    145   MEvtLoop pedloop;
    146   pedloop.SetParList(&plist1);
    147   //pedloop.SetDisplay(display);
    148 
    149   cout << "*************************" << endl;
    150   cout << "** COMPUTING PEDESTALS **" << endl;
    151   cout << "*************************" << endl;
    152 
    153   if (!pedloop.Eventloop())
     122
     123  // If you want to exclude pixels from the beginning, read
     124  // an ascii-file with the corr. pixel numbers (see MBadPixelsCam) 
     125  //badcam.AsciiRead("badpixels.dat");
     126
     127  MJPedestal pedloop;
     128  pedloop.SetInput(&pediter);
     129  pedloop.SetExtractor(&extractor);
     130  pedloop.SetBadPixels(badcam);
     131
     132  if (!pedloop.Process())
    154133    return;
    155  
    156   tlist1.PrintStatistics();
    157 
     134
     135  MPedestalCam pedcam = pedloop.GetPedestalCam();
     136 
    158137  /*****************************/
    159138  /* SECOND LOOP: CALIBRATION  */
    160139  /*****************************/       
    161140
    162   MParList  plist2;
    163   MTaskList tlist2;
    164   plist2.AddToList(&tlist2);
    165   plist2.AddToList(&pedcam);
    166   plist2.AddToList(&badcam);
    167 
    168   // new containers
    169   MGeomCamMagic          geomcam;
    170   MExtractedSignalCam    sigcam;
    171   MCalibrationChargeCam  calcam;
    172   MHCalibrationChargeCam histcharge;
    173   MHCalibrationRelTimeCam     histtime;
    174 
    175 
    176   plist2.AddToList(&geomcam);
    177   plist2.AddToList(&sigcam);
    178   plist2.AddToList(&calcam);
    179   plist2.AddToList(&histcharge);
    180 
    181   //tasks
    182   MReadMarsFile read2("Events");
    183   static_cast<MRead&>(read2).AddFiles(caliter);
    184   read2.DisableAutoScheme();     
    185  
    186   MExtractSignal         sigcalc;
    187   sigcalc.SetRange(hifirst,hilast,lofirst,lolast);
    188   MCalibrationChargeCalc calcalc;
    189   MFCosmics              cosmics;
    190   MContinue              cont(&cosmics);
    191 
    192  
    193   calcalc.SetPulserColor(MCalibrationCam::kCT1);
    194   MFillH fillcam  ("MHCalibrationChargeCam"     , "MExtractedSignalCam");
    195 
    196   tlist2.AddToList(&read2);
    197   tlist2.AddToList(&geomapl);
    198   tlist2.AddToList(&sigcalc);
    199   tlist2.AddToList(&cont);
    200   tlist2.AddToList(&fillcam);
    201   tlist2.AddToList(&calcalc);
    202  
    203   // Create and setup the eventloop
    204   MEvtLoop calloop;
    205   calloop.SetParList(&plist2);
    206   //calloop.SetDisplay(display);
    207 
    208   cout << "***************************" << endl;
    209   cout << "** COMPUTING CALIBRATION **" << endl;
    210   cout << "***************************" << endl;
    211  
    212   if (!calloop.Eventloop())
     141  MJCalibration calloop;
     142
     143  calloop.SetExtractor(&extractor);
     144  calloop.SetInput(&caliter);
     145  calloop.SetQECam(qecam);
     146  calloop.SetBadPixels(pedloop.GetBadPixels());
     147  if (!calloop.Process(pedcam))
    213148    return;
    214  
    215   tlist2.PrintStatistics();
    216  
    217   MLog gauglog;
    218   gauglog.SetOutputFile(Form("%s%s",calcam.GetName(),".txt"),1);
    219   calcam.SetLogStream(&gauglog);
    220   calcam.Print();
    221   calcam.SetLogStream(&gLog);
    222 #if 0
    223   // Create histograms to display
    224   MHCamera disp1  (geomcam, "Cal;Charge",         "Fitted Mean Charges");
    225   MHCamera disp2  (geomcam, "Cal;SigmaCharge",    "Sigma of Fitted Charges");
    226   MHCamera disp3  (geomcam, "Cal;FitProb",        "Probability of Fit");
    227   MHCamera disp4  (geomcam, "Cal;RSigma",         "Reduced Sigmas");
    228   MHCamera disp5  (geomcam, "Cal;RSigma/Charge",  "Reduced Sigma per Charge");
    229   MHCamera disp6  (geomcam, "Cal;FFactorPh",      "Nr. of Photo-electrons (F-Factor Method)");
    230   MHCamera disp7  (geomcam, "Cal;FFactorConv",    "Conversion Factor to photons (F-Factor Method)");
    231   MHCamera disp8  (geomcam, "Cal;FFactorFFactor", "Total F-Factor (F-Factor Method)");
    232   MHCamera disp9  (geomcam, "Cal;BlindPixPh",     "Photon flux inside plexiglass (Blind Pixel Method)");
    233   MHCamera disp10 (geomcam, "Cal;BlindPixConv",   "Conversion Factor to photons (Blind Pixel Method)");
    234   MHCamera disp11 (geomcam, "Cal;BlindPixFFactor","Total F-Factor (Blind Pixel Method)");
    235   MHCamera disp12 (geomcam, "Cal;PINDiodePh",     "Photon flux outside plexiglass (PIN Diode Method)");
    236   MHCamera disp13 (geomcam, "Cal;PINDiodeConv",   "Conversion Factor tp photons (PIN Diode Method)");
    237   MHCamera disp14 (geomcam, "Cal;PINDiodeFFactor","Total F-Factor (PIN Diode Method)");
    238   MHCamera disp15 (geomcam, "Cal;Excluded",       "Pixels previously excluded");
    239   MHCamera disp16 (geomcam, "Cal;NotFitted",      "Pixels that could not be fitted");
    240   MHCamera disp17 (geomcam, "Cal;NotFitValid",    "Pixels with not valid fit results");
    241   MHCamera disp18 (geomcam, "Cal;HiGainOscillating",     "Oscillating Pixels HI Gain");
    242   MHCamera disp19 (geomcam, "Cal;LoGainOscillating",     "Oscillating Pixels LO Gain");
    243   MHCamera disp20 (geomcam, "Cal;HiGainPickup",          "Number Pickup events Hi Gain");
    244   MHCamera disp21 (geomcam, "Cal;LoGainPickup",          "Number Pickup events Lo Gain");
    245   MHCamera disp22 (geomcam, "Cal;Saturation",            "Pixels with saturated Hi Gain");
    246   MHCamera disp23 (geomcam, "Cal;FFactorValid",          "Pixels with valid F-Factor calibration");
    247   MHCamera disp24 (geomcam, "Cal;BlindPixelValid",       "Pixels with valid BlindPixel calibration");
    248   MHCamera disp25 (geomcam, "Cal;PINdiodeFFactorValid",  "Pixels with valid PINDiode calibration");
    249  
    250   MHCamera disp26 (geomcam, "Cal;Ped",         "Pedestals");
    251   MHCamera disp27 (geomcam, "Cal;PedRms",      "Pedestal RMS");
    252  
    253   MHCamera disp28 (geomcam, "time;Time",        "Rel. Arrival Times");
    254   MHCamera disp29 (geomcam, "time;SigmaTime",   "Sigma of Rel. Arrival Times");
    255   MHCamera disp30 (geomcam, "time;TimeProb",    "Probability of Time Fit");
    256   MHCamera disp31 (geomcam, "time;NotFitValid", "Pixels with not valid fit results");
    257   MHCamera disp32 (geomcam, "time;Oscillating", "Oscillating Pixels");
    258  
    259   MHCamera disp33 (geomcam, "Cal;AbsTimeMean", "Abs. Arrival Times");
    260   MHCamera disp34 (geomcam, "Cal;AbsTimeRms",  "RMS of Arrival Times");
    261 
    262   // Fitted charge means and sigmas
    263   disp1.SetCamContent(calcam,  0);
    264   disp1.SetCamError(  calcam,  1);
    265   disp2.SetCamContent(calcam,  2);
    266   disp2.SetCamError(  calcam,  3);
    267 
    268   // Fit probabilities
    269   disp3.SetCamContent(calcam,  4);
    270  
    271   // Reduced Sigmas and reduced sigmas per charge
    272   disp4.SetCamContent(calcam,  5);
    273   disp4.SetCamError(  calcam,  6);
    274   disp5.SetCamContent(calcam,  7);
    275   disp5.SetCamError(  calcam,  8);
    276  
    277   // F-Factor Method
    278   disp6.SetCamContent(calcam,  9);
    279   disp6.SetCamError(  calcam, 10);
    280   disp7.SetCamContent(calcam, 11);
    281   disp7.SetCamError(  calcam, 12);
    282   disp8.SetCamContent(calcam, 13);
    283   disp8.SetCamError(  calcam, 14);
    284  
    285   // Blind Pixel Method
    286   disp9.SetCamContent(calcam, 15);
    287   disp9.SetCamError(  calcam, 16);
    288   disp10.SetCamContent(calcam,17);
    289   disp10.SetCamError(  calcam,18);
    290   disp11.SetCamContent(calcam,19);
    291   disp11.SetCamError(  calcam,20);
    292  
    293   // PIN Diode Method
    294   disp12.SetCamContent(calcam,21);
    295   disp12.SetCamError(  calcam,22);
    296   disp13.SetCamContent(calcam,23);
    297   disp13.SetCamError(  calcam,24);
    298   disp14.SetCamContent(calcam,25);
    299   disp14.SetCamError(  calcam,26);
    300  
    301   // Pixels with defects
    302   disp15.SetCamContent(calcam,27);
    303   disp16.SetCamContent(calcam,28);
    304   disp17.SetCamContent(badcam,9);
    305   disp18.SetCamContent(badcam,16);
    306   disp19.SetCamContent(badcam,15);
    307   disp20.SetCamContent(calcam,29);
    308   disp21.SetCamContent(calcam,30);
    309  
    310   // Lo Gain calibration
    311   disp22.SetCamContent(calcam,31);
    312  
    313   // Valid flags
    314   disp23.SetCamContent(calcam,32);
    315   disp24.SetCamContent(calcam,33);
    316   disp25.SetCamContent(calcam,34);
    317  
    318   // Pedestals
    319   disp26.SetCamContent(calcam,35);
    320   disp26.SetCamError(  calcam,36);
    321   disp27.SetCamContent(calcam,37);
    322   disp27.SetCamError(  calcam,38);
    323  
    324   // Relative Times
    325   disp28.SetCamContent(histtime,0);
    326   disp28.SetCamError(  histtime,1);
    327   disp29.SetCamContent(histtime,2);
    328   disp29.SetCamError(  histtime,3);
    329   disp30.SetCamContent(histtime,4);
    330   disp31.SetCamContent(histtime,5);
    331   disp32.SetCamContent(histtime,6);
    332  
    333   // Absolute Times
    334   disp33.SetCamContent(calcam,39);
    335   disp33.SetCamError(  calcam,40);
    336   disp34.SetCamContent(calcam,41);
    337  
    338   disp1.SetYTitle("Charge [FADC units]");
    339   disp2.SetYTitle("\\sigma_{Charge} [FADC units]");
    340   disp3.SetYTitle("P_{Charge} [1]");
    341  
    342   disp4.SetYTitle("\\sqrt{\\sigma^{2}_{Charge} - RMS^{2}_{Ped}} [FADC Counts]");
    343   disp5.SetYTitle("Reduced Sigma / Mean Charge [1]");
    344  
    345   disp6.SetYTitle("Nr. Photo-electrons [1]");
    346   disp7.SetYTitle("Conversion Factor [Ph/FADC Count]");
    347   disp8.SetYTitle("\\sqrt{N_{Ph}}*\\sigma_{Charge}/\\mu_{Charge} [1] ");
    348  
    349   disp9.SetYTitle("Photon flux [ph/mm^2]");
    350   disp10.SetYTitle("Conversion Factor [Phot/FADC Count]");
    351   disp11.SetYTitle("\\sqrt{N_{Ph}}*\\sigma_{Charge}/\\mu_{Charge} [1]");
    352  
    353   disp12.SetYTitle("Photon flux [ph/mm^2]");
    354   disp13.SetYTitle("Conversion Factor [Phot/FADC Count]");
    355   disp14.SetYTitle("\\sqrt{N_{Ph}}*\\sigma_{Charge}/\\mu_{Charge} [1]");
    356  
    357   disp15.SetYTitle("[1]");
    358   disp16.SetYTitle("[1]");
    359   disp17.SetYTitle("[1]");
    360   disp18.SetYTitle("[1]");
    361   disp19.SetYTitle("[1]");
    362   disp20.SetYTitle("[1]");
    363   disp21.SetYTitle("[1]");
    364   disp22.SetYTitle("[1]");
    365   disp23.SetYTitle("[1]");
    366   disp24.SetYTitle("[1]");
    367   disp25.SetYTitle("[1]");
    368 
    369   disp26.SetYTitle("Ped [FADC Counts ]");
    370   disp27.SetYTitle("RMS_{Ped} [FADC Counts ]");
    371  
    372   disp28.SetYTitle("Time Offset [ns]");
    373   disp29.SetYTitle("Timing resolution [ns]");
    374   disp30.SetYTitle("P_{Time} [1]");
    375  
    376   disp31.SetYTitle("[1]");
    377   disp32.SetYTitle("[1]");
    378  
    379   disp33.SetYTitle("Mean Abs. Time [FADC slice]");
    380   disp34.SetYTitle("RMS Abs. Time [FADC slices]");
    381  
    382   gStyle->SetOptStat(1111);
    383   gStyle->SetOptFit();
    384  
    385   // Charges
    386   TCanvas &c1 = display->AddTab("Fit.Charge");
    387   c1.Divide(2, 3);
    388  
    389   CamDraw(c1, disp1,calcam,1, 2 , 2);
    390   CamDraw(c1, disp2,calcam,2, 2 , 2);
    391 
    392   // Fit Probability
    393   TCanvas &c2 = display->AddTab("Fit.Prob");
    394   c2.Divide(1,3);
    395 
    396   CamDraw(c2, disp3,calcam,1, 1 , 4);
    397  
    398   // Reduced Sigmas
    399   TCanvas &c3 = display->AddTab("Red.Sigma");
    400   c3.Divide(2,3);
    401  
    402   CamDraw(c3, disp4,calcam,1, 2 , 2);
    403   CamDraw(c3, disp5,calcam,2, 2 , 2);
    404 
    405   // F-Factor Method
    406   TCanvas &c4 = display->AddTab("F-Factor");
    407   c4.Divide(3,3);
    408  
    409   CamDraw(c4, disp6,calcam,1, 3 , 2);
    410   CamDraw(c4, disp7,calcam,2, 3 , 2);
    411   CamDraw(c4, disp8,calcam,3, 3 , 2);
    412  
    413   // Blind Pixel Method
    414   TCanvas &c5 = display->AddTab("BlindPix");
    415   c5.Divide(3, 3);
    416  
    417   CamDraw(c5, disp9,calcam,1, 3 ,  9);
    418   CamDraw(c5, disp10,calcam,2, 3 , 2);
    419   CamDraw(c5, disp11,calcam,3, 3 , 2);
    420  
    421   // PIN Diode Method
    422   TCanvas &c6 = display->AddTab("PINDiode");
    423   c6.Divide(3,3);
    424  
    425   CamDraw(c6, disp12,calcam,1, 3 , 9);
    426   CamDraw(c6, disp13,calcam,2, 3 , 2);
    427   CamDraw(c6, disp14,calcam,3, 3 , 2);
    428  
    429   // Defects
    430   TCanvas &c7 = display->AddTab("Defects");
    431   c7.Divide(4,2);
    432  
    433   CamDraw(c7, disp15,calcam,1,4, 0);
    434   CamDraw(c7, disp16,calcam,2,4, 0);
    435   CamDraw(c7, disp20,calcam,3,4, 0);
    436   CamDraw(c7, disp21,calcam,4,4, 0);
    437  
    438   // BadCam
    439   TCanvas &c8 = display->AddTab("Defects");
    440   c8.Divide(3,2);
    441  
    442   CamDraw(c8, disp17,badcam,1,3, 0);
    443   CamDraw(c8, disp18,badcam,2,3, 0);
    444   CamDraw(c8, disp19,badcam,3,3, 0);
    445  
    446   // Valid flags
    447   TCanvas &c9 = display->AddTab("Validity");
    448   c9.Divide(4,2);
    449  
    450   CamDraw(c9, disp22,calcam,1,4,0);
    451   CamDraw(c9, disp23,calcam,2,4,0);
    452   CamDraw(c9, disp24,calcam,3,4,0);
    453   CamDraw(c9, disp25,calcam,4,4,0);
    454  
    455   // Pedestals
    456   TCanvas &c10 = display->AddTab("Pedestals");
    457   c10.Divide(2,3);
    458  
    459   CamDraw(c10,disp26,calcam,1,2,1);
    460   CamDraw(c10,disp27,calcam,2,2,2);
    461  
    462   // Rel. Times
    463   TCanvas &c11 = display->AddTab("Fitted Rel. Times");
    464   c11.Divide(3,3);
    465  
    466   CamDraw(c11,disp28,calcam,1,3,2);
    467   CamDraw(c11,disp29,calcam,2,3,2);
    468   CamDraw(c11,disp30,calcam,3,3,4);
    469  
    470   // Time Defects
    471   TCanvas &c12 = display->AddTab("Time Def.");
    472   c12.Divide(2,2);
    473  
    474   CamDraw(c12, disp31,calcam,1,2, 0);
    475   CamDraw(c12, disp32,calcam,2,2, 0);
    476  
    477   // Abs. Times
    478   TCanvas &c13 = display->AddTab("Abs. Times");
    479   c13.Divide(2,3);
    480  
    481   CamDraw(c13,disp33,calcam,1,2,2);
    482   CamDraw(c13,disp34,calcam,2,2,2);
    483 #endif 
     149
     150  MCalibrationChargeCam calcam = calloop.GetCalibrationCam();
     151
    484152  /************************************************************************/
    485153  /*                THIRD LOOP: PEDESTAL CALIBRATION INTO PHOTONS         */
    … …  
    490158 
    491159  // containers
    492   MCerPhotEvt    nphot;
    493   MPedPhotCam    nphotrms;
     160  MCerPhotEvt         nphot;
     161  MPedPhotCam         nphotrms;
     162  MExtractedSignalCam sigcam;
    494163 
    495164  plist3.AddToList(&geomcam);
    496165  plist3.AddToList(&pedcam);
    497166  plist3.AddToList(&calcam);
    498   //  plist3.AddToList(&badcam);
     167  plist3.AddToList(&badcam);
    499168  plist3.AddToList(&sigcam);
    500169  plist3.AddToList(&nphot);
    … …  
    505174  if(calflag==0)
    506175    calMode=MCalibrate::kNone;
    507  
     176  if(calflag==-1)
     177    calMode=MCalibrate::kDummy;
    508178
    509179  //tasks
    … …  
    511181  static_cast<MRead&>(read3).AddFiles(pediter);
    512182  read3.DisableAutoScheme();
    513   
    514   MExtractSignal  extsig;
    515   extsig.SetRange(hifirst,hilast,lofirst,lolast);
     183 
     184  //  MExtractSignal  extsig;
     185  //  extsig.SetRange(hifirst,hilast,lofirst,lolast);
    516186  MCalibrate      photcalc(calMode);
    517187  MPedPhotCalc    photrmscalc;
    … …  
    519189  tlist3.AddToList(&read3);
    520190  tlist3.AddToList(&geomapl);
    521   tlist3.AddToList(&extsig);
     191  tlist3.AddToList(&extractor);
    522192  tlist3.AddToList(&photcalc);
    523193  tlist3.AddToList(&photrmscalc);
    … …  
    528198  if (!evtloop3.Eventloop())
    529199    return;
    530 
     200 
    531201  tlist3.PrintStatistics();
    532  
     202
    533203  /************************************************************************/
    534204  /*                FOURTH LOOP: DATA CALIBRATION INTO PHOTONS            */
    … …  
    576246  MWriteRootFile write(outname,"RECREATE");
    577247
    578   //  write.AddContainer("MGeomCam"              , "RunHeaders");
    579   //  write.AddContainer("MRawRunHeader"         , "RunHeaders");
    580   //  write.AddContainer("MSrcPosCam"            , "RunHeaders");
    581   //  write.AddContainer("MCalibrationChargeCam" , "RunHeaders");
    582   //  write.AddContainer("MPedPhotCam"           , "RunHeaders"); // Attention, was in Events - Tree!!
    583   //  write.AddContainer("MPedestalCam"          , "RunHeaders");
    584   //  write.AddContainer("MHCalibrationRelTimeCam","RunHeaders");
    585 
    586   //  write.AddContainer("MCerPhotEvt"   , "Events");
    587   //  write.AddContainer("MRawEvtHeader" , "Events");
    588   //  write.AddContainer("MBadPixelsCam" , "Events");
    589   //  write.AddContainer("MPedPhotCam"   , "Events");
    590 
    591   write.AddContainer("MHillas"       , "Parameters");
    592   write.AddContainer("MHillasSrc"    , "Parameters");
    593   write.AddContainer("MHillasExt"    , "Parameters");
    594   write.AddContainer("MNewImagePar"  , "Parameters");
    595   write.AddContainer("MRawEvtHeader" , "Parameters");
    596   write.AddContainer("MRawRunHeader" , "Parameters");
     248  write.AddContainer("MHillas"        , "Parameters");
     249  write.AddContainer("MHillasSrc"     , "Parameters");
     250  write.AddContainer("MHillasExt"     , "Parameters");
     251  write.AddContainer("MNewImagePar"   , "Parameters");
     252  write.AddContainer("MRawEvtHeader"  , "Parameters");
     253  write.AddContainer("MRawRunHeader"  , "Parameters");
    597254  write.AddContainer("MConcentration" , "Parameters");
     255  write.AddContainer("MSrcPosCam"     , "Parameters");
    598256 
    599257  tlist4.AddToList(&read4);
    600258  tlist4.AddToList(&geomapl);
    601   tlist4.AddToList(&extsig);
     259  tlist4.AddToList(&extractor);
    602260  tlist4.AddToList(&photcalc);
    603261  //tlist4.AddToList(&blind);
    … …  
    612270  MEvtLoop datloop;
    613271  datloop.SetParList(&plist4);
    614   //  MProgressBar bar;
    615   //  datloop.SetProgressBar(&bar);
    616272
    617273  cout << "*************************************************************" << endl;
    … …  
    621277  if (!datloop.Eventloop(nmaxevents))
    622278    return;
    623  
     279
    624280  tlist4.PrintStatistics();   
    625   nphotrms.Print();
    626 
    627281}
    628 
    629 #if 0
    630 void CamDraw(TCanvas &c, MHCamera &cam, MCamEvent &evt, Int_t i, Int_t j, Int_t fit)
    631 {
    632 
    633   c.cd(i);
    634   gPad->SetBorderMode(0);
    635   MHCamera *obj1=(MHCamera*)cam.DrawCopy("hist");
    636   //  obj1->AddNotify(evt);
    637  
    638   c.cd(i+j);
    639   gPad->SetBorderMode(0);
    640   obj1->Draw();
    641   ((MHCamera*)obj1)->SetPrettyPalette();
    642 
    643   if (fit != 0)
    644     {
    645       c.cd(i+2*j);
    646       gPad->SetBorderMode(0);
    647       TH1D *obj2 = (TH1D*)obj1->Projection(obj1.GetName());
    648      
    649 //      obj2->Sumw2();
    650       obj2->Draw();
    651       obj2->SetBit(kCanDelete);
    652 
    653       const Double_t min   = obj2->GetBinCenter(obj2->GetXaxis()->GetFirst());
    654       const Double_t max   = obj2->GetBinCenter(obj2->GetXaxis()->GetLast());
    655       const Double_t integ = obj2->Integral("width")/2.5066283;
    656       const Double_t mean  = obj2->GetMean();
    657       const Double_t rms   = obj2->GetRMS();
    658       const Double_t width = max-min;
    659 
    660       if (rms == 0. || width == 0. )
    661         return;
    662      
    663       switch (fit)
    664         {
    665         case 1:
    666           TF1 *sgaus = new TF1("sgaus","gaus(0)",min,max);
    667           sgaus->SetBit(kCanDelete);
    668           sgaus->SetParNames("Area","#mu","#sigma");
    669           sgaus->SetParameters(integ/rms,mean,rms);
    670           sgaus->SetParLimits(0,0.,integ);
    671           sgaus->SetParLimits(1,min,max);
    672           sgaus->SetParLimits(2,0,width/1.5);
    673           obj2->Fit("sgaus","QLR");
    674           obj2->GetFunction("sgaus")->SetLineColor(kYellow);
    675           break;
    676 
    677         case 2:
    678           TString dgausform = "([0]-[3])/[2]*exp(-0.5*(x-[1])*(x-[1])/[2]/[2])";
    679           dgausform += "+[3]/[5]*exp(-0.5*(x-[4])*(x-[4])/[5]/[5])";
    680           TF1 *dgaus = new TF1("dgaus",dgausform.Data(),min,max);
    681           dgaus->SetBit(kCanDelete);
    682           dgaus->SetParNames("A_{tot}","#mu_{1}","#sigma_{1}","A_{2}","#mu_{2}","#sigma_{2}");
    683           dgaus->SetParameters(integ,(min+mean)/2.,width/4.,
    684                                integ/width/2.,(max+mean)/2.,width/4.);
    685           // The left-sided Gauss
    686           dgaus->SetParLimits(0,integ-1.5,integ+1.5);
    687           dgaus->SetParLimits(1,min+(width/10.),mean);
    688           dgaus->SetParLimits(2,0,width/2.);
    689           // The right-sided Gauss
    690           dgaus->SetParLimits(3,0,integ);
    691           dgaus->SetParLimits(4,mean,max-(width/10.));
    692           dgaus->SetParLimits(5,0,width/2.);
    693           obj2->Fit("dgaus","QLRM");
    694           obj2->GetFunction("dgaus")->SetLineColor(kYellow);
    695           break;
    696          
    697         case 3:
    698           TString tgausform = "([0]-[3]-[6])/[2]*exp(-0.5*(x-[1])*(x-[1])/[2]/[2])";
    699           tgausform += "+[3]/[5]*exp(-0.5*(x-[4])*(x-[4])/[5]/[5])";
    700           tgausform += "+[6]/[8]*exp(-0.5*(x-[7])*(x-[7])/[8]/[8])";
    701           TF1 *tgaus = new TF1("tgaus",tgausform.Data(),min,max);
    702           tgaus->SetBit(kCanDelete);
    703           tgaus->SetParNames("A_{tot}","#mu_{1}","#sigma_{1}",
    704                              "A_{2}","#mu_{2}","#sigma_{2}",
    705                              "A_{3}","#mu_{3}","#sigma_{3}");
    706           tgaus->SetParameters(integ,(min+mean)/2,width/4.,
    707                                integ/width/3.,(max+mean)/2.,width/4.,
    708                                integ/width/3.,mean,width/2.);
    709           // The left-sided Gauss
    710           tgaus->SetParLimits(0,integ-1.5,integ+1.5);
    711           tgaus->SetParLimits(1,min+(width/10.),mean);
    712           tgaus->SetParLimits(2,width/15.,width/2.);
    713           // The right-sided Gauss
    714           tgaus->SetParLimits(3,0.,integ);
    715           tgaus->SetParLimits(4,mean,max-(width/10.));
    716           tgaus->SetParLimits(5,width/15.,width/2.);
    717           // The Gauss describing the outliers
    718           tgaus->SetParLimits(6,0.,integ);
    719           tgaus->SetParLimits(7,min,max);
    720           tgaus->SetParLimits(8,width/4.,width/1.5);
    721           obj2->Fit("tgaus","QLRM");
    722           obj2->GetFunction("tgaus")->SetLineColor(kYellow);
    723           break;
    724         case 4:
    725           obj2->Fit("pol0","Q");
    726           obj2->GetFunction("pol0")->SetLineColor(kYellow);
    727           break;
    728         case 9:
    729           break;
    730         default:
    731           obj2->Fit("gaus","Q");
    732           obj2->GetFunction("gaus")->SetLineColor(kYellow);
    733           break;
    734         }
    735      
    736         TArrayI s0(3);
    737         s0[0] = 6;
    738         s0[1] = 1;
    739         s0[2] = 2;
    740 
    741         TArrayI s1(3);
    742         s1[0] = 3;
    743         s1[1] = 4;
    744         s1[2] = 5;
    745 
    746         TArrayI inner(1);
    747         inner[0] = 0;
    748 
    749         TArrayI outer(1);
    750         outer[0] = 1;
    751 
    752         // Just to get the right (maximum) binning
    753         TH1D *half[4];
    754         half[0] = obj1->ProjectionS(s0, inner, "Sector 6-1-2 Inner");
    755         half[1] = obj1->ProjectionS(s1, inner, "Sector 3-4-5 Inner");
    756         half[2] = obj1->ProjectionS(s0, outer, "Sector 6-1-2 Outer");
    757         half[3] = obj1->ProjectionS(s1, outer, "Sector 3-4-5 Outer");
    758 
    759         for (int i=0; i<4; i++)     
    760         {
    761             half[i]->SetLineColor(kRed+i);
    762             half[i]->SetDirectory(0);
    763             half[i]->SetBit(kCanDelete);
    764             half[i]->Draw("same");
    765         }
    766 
    767       gPad->Modified();
    768       gPad->Update();
    769      
    770     }
    771 }
    772 #endif
     282//-------------------------------------------------------------------------------
    773283
    774284Bool_t readDatacards(TString& filename)
Note: See TracChangeset for help on using the changeset viewer.