Changeset 3874 for trunk/MagicSoft


Ignore:
Timestamp:
04/28/04 19:25:36 (20 years ago)
Author:
rico
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars
Files:
7 edited

Legend:

Unmodified
Added
Removed
  • trunk/MagicSoft/Mars/Changelog

    r3873 r3874  
    1818
    1919                                                 -*-*- END OF LINE -*-*-
     20 2004/04/28: Javier Rico
     21
     22   * mbase/MRunIter.[h,cc]
     23     - include function to read runs specified in a character chain
     24       with the format, e.g., run1,run2-run3,run4-run5,...
     25       
     26       
    2027 2004/04/28: Abelardo Moralejo
    2128
  • trunk/MagicSoft/Mars/mbase/MRunIter.cc

    r3002 r3874  
    1717!
    1818!   Author(s): Thomas Bretz, 1/2004 <mailto:tbretz@astro.uni-wuerzburg.de>
     19!   Author(s): Javier Rico,  4/2004 <mailto:jrico@ifae.es>
    1920!
    2021!   Copyright: MAGIC Software Development, 2000-2004
     
    3536
    3637#include <TSystem.h>
     38#include <iostream>
    3739
    3840ClassImp(MRunIter);
     
    6062    return AddDirectory(gSystem->DirName(name), gSystem->BaseName(name), -1);
    6163}
     64
     65// --------------------------------------------------------------------------
     66//
     67// Add runs specified in a character chain with the format:
     68// run1,run2-run3,run4-run5,...
     69// e.g  if runrange="100,105-107,110-112,115"
     70// runs 100,105,106,107,110,111,112 and 115 are included in the iterator list
     71//
     72Int_t MRunIter::AddRuns(const char* runrange, const char* path)
     73{
     74  char* last;
     75  char  chcopy[100];
     76  Int_t lowrun=-1;
     77  Int_t upprun;
     78  Int_t totdir=0;
     79
     80  //  cout << "Analyzing chain " << runrange << " in path " << path << endl;
     81
     82  sprintf(chcopy,"%s",runrange);
     83  last=&chcopy[0];
     84
     85  // loop over the elements of the character chain (break inside the loop)
     86  for(char* cp=last;cp;cp++)
     87    {
     88      // look for a digit, a '-' or a ',' or end of string
     89      char c=*cp;
     90      if(! ((c>='0' && c<='9') || c=='-' || c==',' || c=='\0'))
     91        return totdir;
     92     
     93      // if '-' is found, save the previous number as initial run
     94      if(c=='-' && lowrun<0)
     95        {
     96          char chrun[100];
     97          strncpy(chrun,last,cp-last);
     98          lowrun=atoi(chrun);             
     99          last=cp+1;
     100        }
     101      // if ',' or the end of the string are found, save the previous run or run range
     102      if(c==',' || c=='\0')
     103        {
     104          char chrun[100];
     105          strncpy(chrun,last,cp-last);
     106          upprun=atoi(chrun);
     107          if(lowrun>=0 && lowrun<=upprun)
     108            totdir+=AddRuns(lowrun,upprun,path);
     109          else if(lowrun<0)
     110            totdir+=AddRun(upprun,path);         
     111         
     112          if(c=='\0')
     113            break;
     114          lowrun=-1;
     115          last=cp+1;
     116        }
     117    }
     118  return totdir;
     119}
  • trunk/MagicSoft/Mars/mbase/MRunIter.h

    r2993 r3874  
    4242
    4343    Int_t AddRun(UInt_t run, const char *path=0);
     44    Int_t AddRuns(const char* runrange, const char *path=0);
    4445    Int_t AddRuns(UInt_t from, UInt_t to, const char *path=0)
    4546    {
  • trunk/MagicSoft/Mars/mtemp/mifae/Changelog

    r3859 r3874  
    2020       
    2121 2004/04/28: Javier Rico
     22   * makeHillas.cc, makehillas.datacard
     23     - Use MJPedestal and MJCalibration and simplified datacard
     24
     25   * Makefile
     26     - include mjobs
     27
     28       
     29 2004/04/28: Javier Rico
    2230   * mifae/Changelog
    2331     - Added. Next changes in directory mifae will be reported
  • trunk/MagicSoft/Mars/mtemp/mifae/Makefile

    r3857 r3874  
    3030INCLUDES = -I. \
    3131           -I../../mbase \
     32           -I../../mjobs \
    3233           -I../../mpedestal \
    3334           -I../../mbadpixels \
  • trunk/MagicSoft/Mars/mtemp/mifae/makeHillas.cc

    r3857 r3874  
    1919#include "MHCalibrationChargeCam.h"
    2020#include "MHCalibrationRelTimeCam.h"
    21 #include "MExtractSignal.h"
     21#include "MExtractFixedWindow.h"
    2222#include "MCalibrationChargeCalc.h"
    2323#include "MFCosmics.h"
     
    3939#include "MProgressBar.h"
    4040#include "MArgs.h"
     41#include "MRunIter.h"
     42#include "MJPedestal.h"
     43#include "MJCalibration.h"
    4144
    4245#include <iostream>
    4346#include <fstream>
     47#include <stdlib.h>
    4448
    4549using namespace std;
     
    5660// declaration of variables read from datacards
    5761TString  outname;
    58 TString  pedname;
    59 TString  calname;
    60 TString* datafile;
     62TString  idirname;
     63MRunIter caliter;
     64MRunIter pediter;
     65MRunIter datiter;
    6166ULong_t  nmaxevents=999999999;
    6267Short_t  calflag=1;
     
    6671
    6772const TString defaultcard="input.datacard";
     73
     74/*************************************************************/
     75// makeHillas usage output
    6876static void Usage()
    6977{
     
    7583}
    7684
     85/*************************************************************/
     86// main program
    7787int main(int argc, char **argv)
    7888{
     
    8595    }
    8696
     97  // get name of input datacard file
    8798  TString datacard  = arg.GetArgumentStr(0);
    8899  if(!datacard.Length())
    89100    datacard = defaultcard;
    90101
     102  // read the datacards
    91103  if(!readDatacards(datacard))
    92104    {
     
    94106      return -1;
    95107    }
     108
     109  // make the hillas file
    96110  makeHillas();
    97111}
    98112
     113/*************************************************************/
     114/*************************************************************/
    99115/*************************************************************/
    100116void makeHillas()
    101117{
    102 #if 0
    103   MStatusDisplay *display = new MStatusDisplay;
    104   display->SetUpdateTime(3000);
    105   display->Resize(850,700);
    106  
    107   gStyle->SetOptStat(1111);
    108   gStyle->SetOptFit();
    109 #endif
     118
     119  // set the signal extractor and calibration mode
     120  MExtractFixedWindow extractor;
     121  extractor.SetRange(hifirst,hilast,lofirst,lolast);
     122
     123  MCalibrate::CalibrationMode_t calMode=MCalibrate::kDefault; 
     124  if(calflag==0)
     125    calMode=MCalibrate::kNone;
     126  if(calflag==-1)
     127    calMode=MCalibrate::kDummy;
     128  MCalibrate  calibrate(calMode);
     129 
     130
     131  // general containers
     132  MBadPixelsCam     badcam;
     133  MGeomCamMagic     geomcam;
     134  MGeomApply        geomapl;
     135 
     136  // If you want to exclude pixels from the beginning, read
     137  // an ascii-file with the corr. pixel numbers (see MBadPixelsCam)
     138  //
     139  //  badcam.AsciiRead("badpixels.dat");
     140 
    110141 
    111142  /************************************/
    112143  /* FIRST LOOP: PEDESTAL COMPUTATION */
    113144  /************************************/
    114  
    115   MParList plist1;
    116   MTaskList tlist1;
    117   plist1.AddToList(&tlist1);
    118  
    119   // containers
    120   MPedestalCam   pedcam;
    121   MBadPixelsCam  badcam; 
    122   // badcam.AsciiRead("badpixels.dat");
    123 
    124   plist1.AddToList(&pedcam);
    125   plist1.AddToList(&badcam);
    126 
    127   //tasks
    128   MReadMarsFile read("Events", pedname);
    129   read.DisableAutoScheme();
    130  
    131   MGeomApply     geomapl;
    132   MPedCalcPedRun pedcalc;
    133 
    134   tlist1.AddToList(&read);
    135   tlist1.AddToList(&geomapl);
    136   tlist1.AddToList(&pedcalc);
    137 
    138   // Create and setup the eventloop
    139   MEvtLoop pedloop;
    140   pedloop.SetParList(&plist1);
    141   //pedloop.SetDisplay(display);
     145
     146  MJPedestal pedloop;
     147  pedloop.SetInput(&pediter);
     148  pedloop.SetBadPixels(badcam);
    142149
    143150  cout << "*************************" << endl;
     
    145152  cout << "*************************" << endl;
    146153
    147   if (!pedloop.Eventloop())
     154  if (!pedloop.Process())
    148155    return;
    149  
    150   tlist1.PrintStatistics();
    151156
    152157  /*****************************/
    153158  /* SECOND LOOP: CALIBRATION  */
    154159  /*****************************/       
    155 
    156   MParList  plist2;
    157   MTaskList tlist2;
    158   plist2.AddToList(&tlist2);
    159   plist2.AddToList(&pedcam);
    160   plist2.AddToList(&badcam);
    161 
    162   // new containers
    163   MGeomCamMagic          geomcam;
    164   MExtractedSignalCam    sigcam;
    165   MCalibrationChargeCam  calcam;
    166   MHCalibrationChargeCam histcharge;
    167   MHCalibrationRelTimeCam     histtime;
    168 
    169   plist2.AddToList(&geomcam);
    170   plist2.AddToList(&sigcam);
    171   plist2.AddToList(&calcam);
    172   plist2.AddToList(&histcharge);
    173 
    174   //tasks
    175   MReadMarsFile read2("Events", calname);
    176   read2.DisableAutoScheme();     
    177  
    178   MExtractSignal         sigcalc;
    179   sigcalc.SetRange(hifirst,hilast,lofirst,lolast);
    180   MCalibrationChargeCalc calcalc;
    181   MFCosmics              cosmics;
    182   MContinue              cont(&cosmics);
    183 
    184  
    185   MFillH fillcam  ("MHCalibrationChargeCam"     , "MExtractedSignalCam");
    186 
    187   tlist2.AddToList(&read2);
    188   tlist2.AddToList(&geomapl);
    189   tlist2.AddToList(&sigcalc);
    190   tlist2.AddToList(&cont);
    191   tlist2.AddToList(&fillcam);
    192   tlist2.AddToList(&calcalc);
    193  
    194   // Create and setup the eventloop
    195   MEvtLoop calloop;
    196   calloop.SetParList(&plist2);
    197   //calloop.SetDisplay(display);
    198 
     160  MCalibrationQECam qecam;
     161  MJCalibration     calloop;
     162  calloop.SetInput(&caliter);
     163  calloop.SetExtractor(&extractor);
     164  //
     165  // Set the corr. cams:
     166  //
     167  calloop.SetQECam(qecam);
     168  calloop.SetBadPixels(pedloop.GetBadPixels());
     169 
     170  //
     171  // Apply rel. time calibration: 
     172  //  calloop.SetRelTimeCalibration();
     173 
     174  // Use as arrival time extractor MArrivalTimeCalc2:
     175  //  calloop.SetArrivalTimeLevel(2);
     176 
    199177  cout << "***************************" << endl;
    200178  cout << "** COMPUTING CALIBRATION **" << endl;
    201179  cout << "***************************" << endl;
    202180 
    203   if (!calloop.Eventloop())
     181  if (!calloop.Process(pedloop.GetPedestalCam()))
    204182    return;
    205183 
    206   tlist2.PrintStatistics();
    207  
    208   MLog gauglog;
    209   gauglog.SetOutputFile(Form("%s%s",calcam.GetName(),".txt"),1);
    210   calcam.SetLogStream(&gauglog);
    211   calcam.Print();
    212   calcam.SetLogStream(&gLog);
    213 #if 0
    214   // Create histograms to display
    215   MHCamera disp1  (geomcam, "Cal;Charge",         "Fitted Mean Charges");
    216   MHCamera disp2  (geomcam, "Cal;SigmaCharge",    "Sigma of Fitted Charges");
    217   MHCamera disp3  (geomcam, "Cal;FitProb",        "Probability of Fit");
    218   MHCamera disp4  (geomcam, "Cal;RSigma",         "Reduced Sigmas");
    219   MHCamera disp5  (geomcam, "Cal;RSigma/Charge",  "Reduced Sigma per Charge");
    220   MHCamera disp6  (geomcam, "Cal;FFactorPh",      "Nr. of Photo-electrons (F-Factor Method)");
    221   MHCamera disp7  (geomcam, "Cal;FFactorConv",    "Conversion Factor to photons (F-Factor Method)");
    222   MHCamera disp8  (geomcam, "Cal;FFactorFFactor", "Total F-Factor (F-Factor Method)");
    223   MHCamera disp9  (geomcam, "Cal;BlindPixPh",     "Photon flux inside plexiglass (Blind Pixel Method)");
    224   MHCamera disp10 (geomcam, "Cal;BlindPixConv",   "Conversion Factor to photons (Blind Pixel Method)");
    225   MHCamera disp11 (geomcam, "Cal;BlindPixFFactor","Total F-Factor (Blind Pixel Method)");
    226   MHCamera disp12 (geomcam, "Cal;PINDiodePh",     "Photon flux outside plexiglass (PIN Diode Method)");
    227   MHCamera disp13 (geomcam, "Cal;PINDiodeConv",   "Conversion Factor tp photons (PIN Diode Method)");
    228   MHCamera disp14 (geomcam, "Cal;PINDiodeFFactor","Total F-Factor (PIN Diode Method)");
    229   MHCamera disp15 (geomcam, "Cal;Excluded",       "Pixels previously excluded");
    230   MHCamera disp16 (geomcam, "Cal;NotFitted",      "Pixels that could not be fitted");
    231   MHCamera disp17 (geomcam, "Cal;NotFitValid",    "Pixels with not valid fit results");
    232   MHCamera disp18 (geomcam, "Cal;HiGainOscillating",     "Oscillating Pixels HI Gain");
    233   MHCamera disp19 (geomcam, "Cal;LoGainOscillating",     "Oscillating Pixels LO Gain");
    234   MHCamera disp20 (geomcam, "Cal;HiGainPickup",          "Number Pickup events Hi Gain");
    235   MHCamera disp21 (geomcam, "Cal;LoGainPickup",          "Number Pickup events Lo Gain");
    236   MHCamera disp22 (geomcam, "Cal;Saturation",            "Pixels with saturated Hi Gain");
    237   MHCamera disp23 (geomcam, "Cal;FFactorValid",          "Pixels with valid F-Factor calibration");
    238   MHCamera disp24 (geomcam, "Cal;BlindPixelValid",       "Pixels with valid BlindPixel calibration");
    239   MHCamera disp25 (geomcam, "Cal;PINdiodeFFactorValid",  "Pixels with valid PINDiode calibration");
    240  
    241   MHCamera disp26 (geomcam, "Cal;Ped",         "Pedestals");
    242   MHCamera disp27 (geomcam, "Cal;PedRms",      "Pedestal RMS");
    243  
    244   MHCamera disp28 (geomcam, "time;Time",        "Rel. Arrival Times");
    245   MHCamera disp29 (geomcam, "time;SigmaTime",   "Sigma of Rel. Arrival Times");
    246   MHCamera disp30 (geomcam, "time;TimeProb",    "Probability of Time Fit");
    247   MHCamera disp31 (geomcam, "time;NotFitValid", "Pixels with not valid fit results");
    248   MHCamera disp32 (geomcam, "time;Oscillating", "Oscillating Pixels");
    249  
    250   MHCamera disp33 (geomcam, "Cal;AbsTimeMean", "Abs. Arrival Times");
    251   MHCamera disp34 (geomcam, "Cal;AbsTimeRms",  "RMS of Arrival Times");
    252 
    253   // Fitted charge means and sigmas
    254   disp1.SetCamContent(calcam,  0);
    255   disp1.SetCamError(  calcam,  1);
    256   disp2.SetCamContent(calcam,  2);
    257   disp2.SetCamError(  calcam,  3);
    258 
    259   // Fit probabilities
    260   disp3.SetCamContent(calcam,  4);
    261  
    262   // Reduced Sigmas and reduced sigmas per charge
    263   disp4.SetCamContent(calcam,  5);
    264   disp4.SetCamError(  calcam,  6);
    265   disp5.SetCamContent(calcam,  7);
    266   disp5.SetCamError(  calcam,  8);
    267  
    268   // F-Factor Method
    269   disp6.SetCamContent(calcam,  9);
    270   disp6.SetCamError(  calcam, 10);
    271   disp7.SetCamContent(calcam, 11);
    272   disp7.SetCamError(  calcam, 12);
    273   disp8.SetCamContent(calcam, 13);
    274   disp8.SetCamError(  calcam, 14);
    275  
    276   // Blind Pixel Method
    277   disp9.SetCamContent(calcam, 15);
    278   disp9.SetCamError(  calcam, 16);
    279   disp10.SetCamContent(calcam,17);
    280   disp10.SetCamError(  calcam,18);
    281   disp11.SetCamContent(calcam,19);
    282   disp11.SetCamError(  calcam,20);
    283  
    284   // PIN Diode Method
    285   disp12.SetCamContent(calcam,21);
    286   disp12.SetCamError(  calcam,22);
    287   disp13.SetCamContent(calcam,23);
    288   disp13.SetCamError(  calcam,24);
    289   disp14.SetCamContent(calcam,25);
    290   disp14.SetCamError(  calcam,26);
    291  
    292   // Pixels with defects
    293   disp15.SetCamContent(calcam,27);
    294   disp16.SetCamContent(calcam,28);
    295   disp17.SetCamContent(badcam,9);
    296   disp18.SetCamContent(badcam,16);
    297   disp19.SetCamContent(badcam,15);
    298   disp20.SetCamContent(calcam,29);
    299   disp21.SetCamContent(calcam,30);
    300  
    301   // Lo Gain calibration
    302   disp22.SetCamContent(calcam,31);
    303  
    304   // Valid flags
    305   disp23.SetCamContent(calcam,32);
    306   disp24.SetCamContent(calcam,33);
    307   disp25.SetCamContent(calcam,34);
    308  
    309   // Pedestals
    310   disp26.SetCamContent(calcam,35);
    311   disp26.SetCamError(  calcam,36);
    312   disp27.SetCamContent(calcam,37);
    313   disp27.SetCamError(  calcam,38);
    314  
    315   // Relative Times
    316   disp28.SetCamContent(histtime,0);
    317   disp28.SetCamError(  histtime,1);
    318   disp29.SetCamContent(histtime,2);
    319   disp29.SetCamError(  histtime,3);
    320   disp30.SetCamContent(histtime,4);
    321   disp31.SetCamContent(histtime,5);
    322   disp32.SetCamContent(histtime,6);
    323  
    324   // Absolute Times
    325   disp33.SetCamContent(calcam,39);
    326   disp33.SetCamError(  calcam,40);
    327   disp34.SetCamContent(calcam,41);
    328  
    329   disp1.SetYTitle("Charge [FADC units]");
    330   disp2.SetYTitle("\\sigma_{Charge} [FADC units]");
    331   disp3.SetYTitle("P_{Charge} [1]");
    332  
    333   disp4.SetYTitle("\\sqrt{\\sigma^{2}_{Charge} - RMS^{2}_{Ped}} [FADC Counts]");
    334   disp5.SetYTitle("Reduced Sigma / Mean Charge [1]");
    335  
    336   disp6.SetYTitle("Nr. Photo-electrons [1]");
    337   disp7.SetYTitle("Conversion Factor [Ph/FADC Count]");
    338   disp8.SetYTitle("\\sqrt{N_{Ph}}*\\sigma_{Charge}/\\mu_{Charge} [1] ");
    339  
    340   disp9.SetYTitle("Photon flux [ph/mm^2]");
    341   disp10.SetYTitle("Conversion Factor [Phot/FADC Count]");
    342   disp11.SetYTitle("\\sqrt{N_{Ph}}*\\sigma_{Charge}/\\mu_{Charge} [1]");
    343  
    344   disp12.SetYTitle("Photon flux [ph/mm^2]");
    345   disp13.SetYTitle("Conversion Factor [Phot/FADC Count]");
    346   disp14.SetYTitle("\\sqrt{N_{Ph}}*\\sigma_{Charge}/\\mu_{Charge} [1]");
    347  
    348   disp15.SetYTitle("[1]");
    349   disp16.SetYTitle("[1]");
    350   disp17.SetYTitle("[1]");
    351   disp18.SetYTitle("[1]");
    352   disp19.SetYTitle("[1]");
    353   disp20.SetYTitle("[1]");
    354   disp21.SetYTitle("[1]");
    355   disp22.SetYTitle("[1]");
    356   disp23.SetYTitle("[1]");
    357   disp24.SetYTitle("[1]");
    358   disp25.SetYTitle("[1]");
    359 
    360   disp26.SetYTitle("Ped [FADC Counts ]");
    361   disp27.SetYTitle("RMS_{Ped} [FADC Counts ]");
    362  
    363   disp28.SetYTitle("Time Offset [ns]");
    364   disp29.SetYTitle("Timing resolution [ns]");
    365   disp30.SetYTitle("P_{Time} [1]");
    366  
    367   disp31.SetYTitle("[1]");
    368   disp32.SetYTitle("[1]");
    369  
    370   disp33.SetYTitle("Mean Abs. Time [FADC slice]");
    371   disp34.SetYTitle("RMS Abs. Time [FADC slices]");
    372  
    373   gStyle->SetOptStat(1111);
    374   gStyle->SetOptFit();
    375  
    376   // Charges
    377   TCanvas &c1 = display->AddTab("Fit.Charge");
    378   c1.Divide(2, 3);
    379  
    380   CamDraw(c1, disp1,calcam,1, 2 , 2);
    381   CamDraw(c1, disp2,calcam,2, 2 , 2);
    382 
    383   // Fit Probability
    384   TCanvas &c2 = display->AddTab("Fit.Prob");
    385   c2.Divide(1,3);
    386 
    387   CamDraw(c2, disp3,calcam,1, 1 , 4);
    388  
    389   // Reduced Sigmas
    390   TCanvas &c3 = display->AddTab("Red.Sigma");
    391   c3.Divide(2,3);
    392  
    393   CamDraw(c3, disp4,calcam,1, 2 , 2);
    394   CamDraw(c3, disp5,calcam,2, 2 , 2);
    395 
    396   // F-Factor Method
    397   TCanvas &c4 = display->AddTab("F-Factor");
    398   c4.Divide(3,3);
    399  
    400   CamDraw(c4, disp6,calcam,1, 3 , 2);
    401   CamDraw(c4, disp7,calcam,2, 3 , 2);
    402   CamDraw(c4, disp8,calcam,3, 3 , 2);
    403  
    404   // Blind Pixel Method
    405   TCanvas &c5 = display->AddTab("BlindPix");
    406   c5.Divide(3, 3);
    407  
    408   CamDraw(c5, disp9,calcam,1, 3 ,  9);
    409   CamDraw(c5, disp10,calcam,2, 3 , 2);
    410   CamDraw(c5, disp11,calcam,3, 3 , 2);
    411  
    412   // PIN Diode Method
    413   TCanvas &c6 = display->AddTab("PINDiode");
    414   c6.Divide(3,3);
    415  
    416   CamDraw(c6, disp12,calcam,1, 3 , 9);
    417   CamDraw(c6, disp13,calcam,2, 3 , 2);
    418   CamDraw(c6, disp14,calcam,3, 3 , 2);
    419  
    420   // Defects
    421   TCanvas &c7 = display->AddTab("Defects");
    422   c7.Divide(4,2);
    423  
    424   CamDraw(c7, disp15,calcam,1,4, 0);
    425   CamDraw(c7, disp16,calcam,2,4, 0);
    426   CamDraw(c7, disp20,calcam,3,4, 0);
    427   CamDraw(c7, disp21,calcam,4,4, 0);
    428  
    429   // BadCam
    430   TCanvas &c8 = display->AddTab("Defects");
    431   c8.Divide(3,2);
    432  
    433   CamDraw(c8, disp17,badcam,1,3, 0);
    434   CamDraw(c8, disp18,badcam,2,3, 0);
    435   CamDraw(c8, disp19,badcam,3,3, 0);
    436  
    437   // Valid flags
    438   TCanvas &c9 = display->AddTab("Validity");
    439   c9.Divide(4,2);
    440  
    441   CamDraw(c9, disp22,calcam,1,4,0);
    442   CamDraw(c9, disp23,calcam,2,4,0);
    443   CamDraw(c9, disp24,calcam,3,4,0);
    444   CamDraw(c9, disp25,calcam,4,4,0);
    445  
    446   // Pedestals
    447   TCanvas &c10 = display->AddTab("Pedestals");
    448   c10.Divide(2,3);
    449  
    450   CamDraw(c10,disp26,calcam,1,2,1);
    451   CamDraw(c10,disp27,calcam,2,2,2);
    452  
    453   // Rel. Times
    454   TCanvas &c11 = display->AddTab("Fitted Rel. Times");
    455   c11.Divide(3,3);
    456  
    457   CamDraw(c11,disp28,calcam,1,3,2);
    458   CamDraw(c11,disp29,calcam,2,3,2);
    459   CamDraw(c11,disp30,calcam,3,3,4);
    460  
    461   // Time Defects
    462   TCanvas &c12 = display->AddTab("Time Def.");
    463   c12.Divide(2,2);
    464  
    465   CamDraw(c12, disp31,calcam,1,2, 0);
    466   CamDraw(c12, disp32,calcam,2,2, 0);
    467  
    468   // Abs. Times
    469   TCanvas &c13 = display->AddTab("Abs. Times");
    470   c13.Divide(2,3);
    471  
    472   CamDraw(c13,disp33,calcam,1,2,2);
    473   CamDraw(c13,disp34,calcam,2,2,2);
    474 #endif 
     184
    475185  /************************************************************************/
    476186  /*                THIRD LOOP: PEDESTAL CALIBRATION INTO PHOTONS         */
     
    485195 
    486196  plist3.AddToList(&geomcam);
    487   plist3.AddToList(&pedcam);
    488   plist3.AddToList(&calcam);
    489   plist3.AddToList(&badcam);
    490   plist3.AddToList(&sigcam);
     197  plist3.AddToList(&pedloop.GetPedestalCam());
     198  plist3.AddToList(&calloop.GetCalibrationCam());
     199  plist3.AddToList(&calloop.GetQECam());
     200  plist3.AddToList(&calloop.GetRelTimeCam());
     201  plist3.AddToList(&calloop.GetBadPixels());
    491202  plist3.AddToList(&nphot);
    492203  plist3.AddToList(&nphotrms);
    493204
    494  
    495   MCalibrate::CalibrationMode_t calMode=MCalibrate::kDefault; 
    496   if(calflag==0)
    497     calMode=MCalibrate::kNone;
    498  
    499 
    500205  //tasks
    501   MReadMarsFile read3("Events", pedname);
     206  MReadMarsFile read3("Events");
    502207  read3.DisableAutoScheme();
    503  
    504   MExtractSignal  extsig;
    505   extsig.SetRange(hifirst,hilast,lofirst,lolast);
    506   MCalibrate      photcalc(calMode);
     208  static_cast<MRead&>(read3).AddFiles(pediter);
     209 
    507210  MPedPhotCalc    photrmscalc;
    508211 
    509212  tlist3.AddToList(&read3);
    510213  tlist3.AddToList(&geomapl);
    511   tlist3.AddToList(&extsig);
    512   tlist3.AddToList(&photcalc);
     214  tlist3.AddToList(&extractor);
     215  tlist3.AddToList(&calibrate);
    513216  tlist3.AddToList(&photrmscalc);
    514217
     
    531234  // containers
    532235  MHillas       hillas;
    533   MSrcPosCam    source;
     236  MSrcPosCam    srcposcam;
    534237  MRawRunHeader runhead;
    535  
     238
    536239  plist4.AddToList(&geomcam);
    537   plist4.AddToList(&pedcam);
    538   plist4.AddToList(&calcam);
    539   plist4.AddToList(&badcam);
     240  plist4.AddToList(&pedloop.GetPedestalCam());
     241  plist4.AddToList(&calloop.GetCalibrationCam());
     242  plist4.AddToList(&calloop.GetQECam());
     243  plist4.AddToList(&calloop.GetRelTimeCam());
     244  plist4.AddToList(&calloop.GetBadPixels());
    540245  plist4.AddToList(&nphot);
    541246  plist4.AddToList(&nphotrms);
    542   plist4.AddToList(&source);
     247  plist4.AddToList(&srcposcam);
    543248  plist4.AddToList(&hillas);
    544249  plist4.AddToList(&runhead);
     
    546251  //tasks
    547252  MReadMarsFile read4("Events");
    548  
    549   for(Int_t i=0;i<nfiles;i++)
    550     read4.AddFile(datafile[i]);
    551253  read4.DisableAutoScheme();
    552  
    553   // set bad pixels
    554   MBlindPixelCalc   blind;
    555   MBlindPixelCalc   blind2;
    556   const Short_t x[16] = {330,395,329,396,389,
    557                          323,388,322,384,385,
    558                          386,387,321,320,319,
    559                          394};
    560   const TArrayS bp(16,(Short_t*)x);
    561   blind.SetPixelIndices(bp);
    562   blind2.SetPixelIndices(bp);
     254  static_cast<MRead&>(read4).AddFiles(datiter); 
    563255 
    564256  MImgCleanStd      clean(lcore,ltail);
     
    568260  MWriteRootFile write(outname,"RECREATE");
    569261
    570   //  write.AddContainer("MGeomCam"              , "RunHeaders");
    571   //  write.AddContainer("MRawRunHeader"         , "RunHeaders");
    572   //  write.AddContainer("MSrcPosCam"            , "RunHeaders");
    573   //  write.AddContainer("MCalibrationChargeCam" , "RunHeaders");
    574   //  write.AddContainer("MPedPhotCam"           , "RunHeaders"); // Attention, was in Events - Tree!!
    575   //  write.AddContainer("MPedestalCam"          , "RunHeaders");
    576   //  write.AddContainer("MHCalibrationRelTimeCam","RunHeaders");
    577 
    578   //  write.AddContainer("MCerPhotEvt"   , "Events");
    579   //  write.AddContainer("MRawEvtHeader" , "Events");
    580   //  write.AddContainer("MBadPixelsCam" , "Events");
    581   //  write.AddContainer("MPedPhotCam"   , "Events");
    582 
    583   write.AddContainer("MHillas"       , "Parameters");
    584   write.AddContainer("MHillasSrc"    , "Parameters");
    585   write.AddContainer("MHillasExt"    , "Parameters");
    586   write.AddContainer("MNewImagePar"  , "Parameters");
    587   write.AddContainer("MRawEvtHeader" , "Parameters");
    588   write.AddContainer("MRawRunHeader" , "Parameters");
     262  write.AddContainer("MHillas"       ,  "Parameters");
     263  write.AddContainer("MHillasSrc"    ,  "Parameters");
     264  write.AddContainer("MHillasExt"    ,  "Parameters");
     265  write.AddContainer("MNewImagePar"  ,  "Parameters");
     266  write.AddContainer("MRawEvtHeader" ,  "Parameters");
     267  write.AddContainer("MRawRunHeader" ,  "Parameters");
    589268  write.AddContainer("MConcentration" , "Parameters");
     269  write.AddContainer("MSrcPosCam",      "Parameters");
    590270 
    591271  tlist4.AddToList(&read4);
    592272  tlist4.AddToList(&geomapl);
    593   tlist4.AddToList(&extsig);
    594   tlist4.AddToList(&photcalc);
    595   //tlist4.AddToList(&blind);
     273  tlist4.AddToList(&extractor);
     274  tlist4.AddToList(&calibrate);
    596275  tlist4.AddToList(&clean);
    597   //tlist4.AddToList(&blind2);
    598276  tlist4.AddToList(&hcalc);
    599   //  tlist4.AddToList(&srcposcalc);
    600277  tlist4.AddToList(&csrc1);
    601278  tlist4.AddToList(&write);
     
    604281  MEvtLoop datloop;
    605282  datloop.SetParList(&plist4);
    606   //  MProgressBar bar;
    607   //  datloop.SetProgressBar(&bar);
    608283
    609284  cout << "*************************************************************" << endl;
     
    615290 
    616291  tlist4.PrintStatistics();   
    617 
    618292}
    619293
    620 #if 0
    621 void CamDraw(TCanvas &c, MHCamera &cam, MCamEvent &evt, Int_t i, Int_t j, Int_t fit)
    622 {
    623 
    624   c.cd(i);
    625   gPad->SetBorderMode(0);
    626   MHCamera *obj1=(MHCamera*)cam.DrawCopy("hist");
    627   //  obj1->AddNotify(evt);
    628  
    629   c.cd(i+j);
    630   gPad->SetBorderMode(0);
    631   obj1->Draw();
    632   ((MHCamera*)obj1)->SetPrettyPalette();
    633 
    634   if (fit != 0)
    635     {
    636       c.cd(i+2*j);
    637       gPad->SetBorderMode(0);
    638       TH1D *obj2 = (TH1D*)obj1->Projection(obj1.GetName());
    639      
    640 //      obj2->Sumw2();
    641       obj2->Draw();
    642       obj2->SetBit(kCanDelete);
    643 
    644       const Double_t min   = obj2->GetBinCenter(obj2->GetXaxis()->GetFirst());
    645       const Double_t max   = obj2->GetBinCenter(obj2->GetXaxis()->GetLast());
    646       const Double_t integ = obj2->Integral("width")/2.5066283;
    647       const Double_t mean  = obj2->GetMean();
    648       const Double_t rms   = obj2->GetRMS();
    649       const Double_t width = max-min;
    650 
    651       if (rms == 0. || width == 0. )
    652         return;
    653      
    654       switch (fit)
    655         {
    656         case 1:
    657           TF1 *sgaus = new TF1("sgaus","gaus(0)",min,max);
    658           sgaus->SetBit(kCanDelete);
    659           sgaus->SetParNames("Area","#mu","#sigma");
    660           sgaus->SetParameters(integ/rms,mean,rms);
    661           sgaus->SetParLimits(0,0.,integ);
    662           sgaus->SetParLimits(1,min,max);
    663           sgaus->SetParLimits(2,0,width/1.5);
    664           obj2->Fit("sgaus","QLR");
    665           obj2->GetFunction("sgaus")->SetLineColor(kYellow);
    666           break;
    667 
    668         case 2:
    669           TString dgausform = "([0]-[3])/[2]*exp(-0.5*(x-[1])*(x-[1])/[2]/[2])";
    670           dgausform += "+[3]/[5]*exp(-0.5*(x-[4])*(x-[4])/[5]/[5])";
    671           TF1 *dgaus = new TF1("dgaus",dgausform.Data(),min,max);
    672           dgaus->SetBit(kCanDelete);
    673           dgaus->SetParNames("A_{tot}","#mu_{1}","#sigma_{1}","A_{2}","#mu_{2}","#sigma_{2}");
    674           dgaus->SetParameters(integ,(min+mean)/2.,width/4.,
    675                                integ/width/2.,(max+mean)/2.,width/4.);
    676           // The left-sided Gauss
    677           dgaus->SetParLimits(0,integ-1.5,integ+1.5);
    678           dgaus->SetParLimits(1,min+(width/10.),mean);
    679           dgaus->SetParLimits(2,0,width/2.);
    680           // The right-sided Gauss
    681           dgaus->SetParLimits(3,0,integ);
    682           dgaus->SetParLimits(4,mean,max-(width/10.));
    683           dgaus->SetParLimits(5,0,width/2.);
    684           obj2->Fit("dgaus","QLRM");
    685           obj2->GetFunction("dgaus")->SetLineColor(kYellow);
    686           break;
    687          
    688         case 3:
    689           TString tgausform = "([0]-[3]-[6])/[2]*exp(-0.5*(x-[1])*(x-[1])/[2]/[2])";
    690           tgausform += "+[3]/[5]*exp(-0.5*(x-[4])*(x-[4])/[5]/[5])";
    691           tgausform += "+[6]/[8]*exp(-0.5*(x-[7])*(x-[7])/[8]/[8])";
    692           TF1 *tgaus = new TF1("tgaus",tgausform.Data(),min,max);
    693           tgaus->SetBit(kCanDelete);
    694           tgaus->SetParNames("A_{tot}","#mu_{1}","#sigma_{1}",
    695                              "A_{2}","#mu_{2}","#sigma_{2}",
    696                              "A_{3}","#mu_{3}","#sigma_{3}");
    697           tgaus->SetParameters(integ,(min+mean)/2,width/4.,
    698                                integ/width/3.,(max+mean)/2.,width/4.,
    699                                integ/width/3.,mean,width/2.);
    700           // The left-sided Gauss
    701           tgaus->SetParLimits(0,integ-1.5,integ+1.5);
    702           tgaus->SetParLimits(1,min+(width/10.),mean);
    703           tgaus->SetParLimits(2,width/15.,width/2.);
    704           // The right-sided Gauss
    705           tgaus->SetParLimits(3,0.,integ);
    706           tgaus->SetParLimits(4,mean,max-(width/10.));
    707           tgaus->SetParLimits(5,width/15.,width/2.);
    708           // The Gauss describing the outliers
    709           tgaus->SetParLimits(6,0.,integ);
    710           tgaus->SetParLimits(7,min,max);
    711           tgaus->SetParLimits(8,width/4.,width/1.5);
    712           obj2->Fit("tgaus","QLRM");
    713           obj2->GetFunction("tgaus")->SetLineColor(kYellow);
    714           break;
    715         case 4:
    716           obj2->Fit("pol0","Q");
    717           obj2->GetFunction("pol0")->SetLineColor(kYellow);
    718           break;
    719         case 9:
    720           break;
    721         default:
    722           obj2->Fit("gaus","Q");
    723           obj2->GetFunction("gaus")->SetLineColor(kYellow);
    724           break;
    725         }
    726      
    727         TArrayI s0(3);
    728         s0[0] = 6;
    729         s0[1] = 1;
    730         s0[2] = 2;
    731 
    732         TArrayI s1(3);
    733         s1[0] = 3;
    734         s1[1] = 4;
    735         s1[2] = 5;
    736 
    737         TArrayI inner(1);
    738         inner[0] = 0;
    739 
    740         TArrayI outer(1);
    741         outer[0] = 1;
    742 
    743         // Just to get the right (maximum) binning
    744         TH1D *half[4];
    745         half[0] = obj1->ProjectionS(s0, inner, "Sector 6-1-2 Inner");
    746         half[1] = obj1->ProjectionS(s1, inner, "Sector 3-4-5 Inner");
    747         half[2] = obj1->ProjectionS(s0, outer, "Sector 6-1-2 Outer");
    748         half[3] = obj1->ProjectionS(s1, outer, "Sector 3-4-5 Outer");
    749 
    750         for (int i=0; i<4; i++)     
    751         {
    752             half[i]->SetLineColor(kRed+i);
    753             half[i]->SetDirectory(0);
    754             half[i]->SetBit(kCanDelete);
    755             half[i]->Draw("same");
    756         }
    757 
    758       gPad->Modified();
    759       gPad->Update();
    760      
    761     }
    762 }
    763 #endif
    764 
     294/******************************************************************************/
    765295Bool_t readDatacards(TString& filename)
    766296{
     
    773303
    774304  TString word;
    775   Int_t ifile = 0;
    776305 
    777306  while(ifun >> word)
     
    788317        ifun >> nmaxevents;
    789318
    790       // pedestal file name
    791       if(strcmp(word.Data(),"PFILE")==0)
    792         {
    793           if(pedname.Length())
    794             cout << "readDataCards Warning: overriding pedestal file name" << endl;
    795           ifun >> pedname;
    796         }
    797 
    798       // calibration file name
    799       if(strcmp(word.Data(),"CFILE")==0)
    800         {
    801           if(calname.Length())
    802             cout << "readDataCards Warning: overriding calibration file name" << endl;
    803           ifun >> calname;
    804         }
    805 
    806       // number of data files
    807       if(strcmp(word.Data(),"NFILES")==0)
    808         {
    809           if(nfiles)
    810             {
    811               Int_t dummy;
    812               cout << "readDataCards Warning: trying to set a new value of number of files. Skiping..." << endl;
    813               ifun >> dummy;
    814             }
    815           else
    816             {
    817               ifun >> nfiles;
    818               datafile = new TString[nfiles];
    819             }
    820         }
    821 
    822 
    823       // data file name
    824       if(strcmp(word.Data(),"DFILE")==0)
    825         {
    826           TString newfile;
    827           ifun >> newfile;
    828           if(ifile<nfiles)
    829             datafile[ifile++]=newfile;
    830           else
    831             {
    832               cout << "readDataCards Error: trying to add more data files than specified" << endl;
    833               return kFALSE;
    834             }
     319
     320      // input file directory
     321      if(strcmp(word.Data(),"IDIR")==0)
     322        {
     323          if(idirname.Length())
     324            cout << "readDataCards Warning: overriding input directory file name" << endl;
     325          ifun >> idirname;
     326        }
     327
     328      // pedestal runs
     329      if(strcmp(word.Data(),"PRUNS")==0)
     330        {
     331          if(pediter.GetNumRuns())
     332            cout << "readDataCards Warning: adding pedestal runs to the existing list" << endl;
     333          ifun >> word;
     334          pediter.AddRuns(word.Data(),idirname.Data());
     335        }
     336
     337      // calibration runs
     338      if(strcmp(word.Data(),"CRUNS")==0)
     339        {
     340          if(caliter.GetNumRuns())
     341            cout << "readDataCards Warning: adding calibration runs to the existing list" << endl;
     342          ifun >> word;
     343          caliter.AddRuns(word.Data(),idirname.Data());
     344        }
     345
     346      // data runs
     347      if(strcmp(word.Data(),"DRUNS")==0)
     348        {
     349          if(datiter.GetNumRuns())
     350            cout << "readDataCards Warning: adding data runs to the existing list" << endl;
     351          ifun >> word;
     352          datiter.AddRuns(word.Data(),idirname.Data());
    835353        }
    836354     
     
    855373    }
    856374
     375  pediter.Reset();
     376  caliter.Reset();
     377  datiter.Reset();
     378  TString pfile;
    857379
    858380  // Dump read values
     
    860382  cout << "* Datacards read from file " << filename << endl;
    861383  cout << "************************************************" << endl;
    862   cout << "Pedestal file: " << pedname << endl;
    863   cout << "Calibration file: " << calname << endl;
    864   cout << "Number of data files: " << nfiles << endl;
    865   cout << "Data files: " << endl;
    866   for(int i=0;i<nfiles;i++)
    867     cout << datafile[i] << endl;
     384  cout << "Pedestal file (s): "  << endl;
     385  while(!(pfile=pediter.Next()).IsNull())
     386    cout << pfile << endl;
     387  cout << "Calibration file (s): "  << endl;
     388  while(!(pfile=caliter.Next()).IsNull())
     389    cout << pfile << endl;
     390  cout << "Data file (s): "  << endl;
     391  while(!(pfile=datiter.Next()).IsNull())
     392    cout << pfile << endl;
    868393  cout << "Maximum number of events: " << nmaxevents << endl;
    869394  cout << "Output file name: " << outname << endl;
     
    872397  cout << "***********" << endl << endl;
    873398
    874   if(!pedname.Length())
     399  if(!pediter.GetNumEntries())
    875400    {
    876401      cout << "No pedestal file name specified" << endl;
    877402      return kFALSE;
    878403    }
    879   if(!calname.Length())
     404  if(!caliter.GetNumEntries())
    880405    {
    881406      cout << "No calibration file name specified" << endl;
  • trunk/MagicSoft/Mars/mtemp/mifae/makehillas.datacard

    r3857 r3874  
    33NEVENTS 99999999
    44
    5 // pedestal and calibration files
    6 PFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16743_P_CrabOn_E.root
    7 CFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16744_C_CrabOn_E.root
     5// data file directory
     6IDIR /local_disk/jrico/rootdata/Crab20040215/
    87
    9 // data files
    10 NFILES 23
    11 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16745_D_CrabOn_E.root
    12 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16746_D_CrabOn_E.root
    13 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16747_D_CrabOn_E.root
    14 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16748_D_CrabOn_E.root
    15 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16749_D_CrabOn_E.root
    16 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16750_D_CrabOn_E.root
    17 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16751_D_CrabOn_E.root
    18 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16752_D_CrabOn_E.root
    19 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16753_D_CrabOn_E.root
    20 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16754_D_CrabOn_E.root
    21 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16755_D_CrabOn_E.root
    22 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16756_D_CrabOn_E.root
    23 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16757_D_CrabOn_E.root
    24 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16758_D_CrabOn_E.root
    25 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16759_D_CrabOn_E.root
    26 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16760_D_CrabOn_E.root
    27 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16761_D_CrabOn_E.root
    28 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16762_D_CrabOn_E.root
    29 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16763_D_CrabOn_E.root
    30 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16764_D_CrabOn_E.root
    31 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16765_D_CrabOn_E.root
    32 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16766_D_CrabOn_E.root
    33 DFILE /local_disk/jrico/rootdata/Crab20040215/20040215_16767_D_CrabOn_E.root
     8// Pedestal (PRUNS), calibration (CRUNS) and data runs (DRUNS), e.g 1500-23444,25444,25455-26544
     9PRUNS 16743
     10CRUNS 16744
     11DRUNS 16745-16767
    3412
    3513// output file name
Note: See TracChangeset for help on using the changeset viewer.