Index: trunk/MagicSoft/Mars/Changelog =================================================================== --- trunk/MagicSoft/Mars/Changelog (revision 3873) +++ trunk/MagicSoft/Mars/Changelog (revision 3874) @@ -18,4 +18,11 @@ -*-*- END OF LINE -*-*- + 2004/04/28: Javier Rico + + * mbase/MRunIter.[h,cc] + - include function to read runs specified in a character chain + with the format, e.g., run1,run2-run3,run4-run5,... + + 2004/04/28: Abelardo Moralejo Index: trunk/MagicSoft/Mars/mbase/MRunIter.cc =================================================================== --- trunk/MagicSoft/Mars/mbase/MRunIter.cc (revision 3873) +++ trunk/MagicSoft/Mars/mbase/MRunIter.cc (revision 3874) @@ -17,4 +17,5 @@ ! ! Author(s): Thomas Bretz, 1/2004 +! Author(s): Javier Rico, 4/2004 ! ! Copyright: MAGIC Software Development, 2000-2004 @@ -35,4 +36,5 @@ #include +#include ClassImp(MRunIter); @@ -60,2 +62,58 @@ return AddDirectory(gSystem->DirName(name), gSystem->BaseName(name), -1); } + +// -------------------------------------------------------------------------- +// +// Add runs specified in a character chain with the format: +// run1,run2-run3,run4-run5,... +// e.g if runrange="100,105-107,110-112,115" +// runs 100,105,106,107,110,111,112 and 115 are included in the iterator list +// +Int_t MRunIter::AddRuns(const char* runrange, const char* path) +{ + char* last; + char chcopy[100]; + Int_t lowrun=-1; + Int_t upprun; + Int_t totdir=0; + + // cout << "Analyzing chain " << runrange << " in path " << path << endl; + + sprintf(chcopy,"%s",runrange); + last=&chcopy[0]; + + // loop over the elements of the character chain (break inside the loop) + for(char* cp=last;cp;cp++) + { + // look for a digit, a '-' or a ',' or end of string + char c=*cp; + if(! ((c>='0' && c<='9') || c=='-' || c==',' || c=='\0')) + return totdir; + + // if '-' is found, save the previous number as initial run + if(c=='-' && lowrun<0) + { + char chrun[100]; + strncpy(chrun,last,cp-last); + lowrun=atoi(chrun); + last=cp+1; + } + // if ',' or the end of the string are found, save the previous run or run range + if(c==',' || c=='\0') + { + char chrun[100]; + strncpy(chrun,last,cp-last); + upprun=atoi(chrun); + if(lowrun>=0 && lowrun<=upprun) + totdir+=AddRuns(lowrun,upprun,path); + else if(lowrun<0) + totdir+=AddRun(upprun,path); + + if(c=='\0') + break; + lowrun=-1; + last=cp+1; + } + } + return totdir; +} Index: trunk/MagicSoft/Mars/mbase/MRunIter.h =================================================================== --- trunk/MagicSoft/Mars/mbase/MRunIter.h (revision 3873) +++ trunk/MagicSoft/Mars/mbase/MRunIter.h (revision 3874) @@ -42,4 +42,5 @@ Int_t AddRun(UInt_t run, const char *path=0); + Int_t AddRuns(const char* runrange, const char *path=0); Int_t AddRuns(UInt_t from, UInt_t to, const char *path=0) { Index: trunk/MagicSoft/Mars/mtemp/mifae/Changelog =================================================================== --- trunk/MagicSoft/Mars/mtemp/mifae/Changelog (revision 3873) +++ trunk/MagicSoft/Mars/mtemp/mifae/Changelog (revision 3874) @@ -20,4 +20,12 @@ 2004/04/28: Javier Rico + * makeHillas.cc, makehillas.datacard + - Use MJPedestal and MJCalibration and simplified datacard + + * Makefile + - include mjobs + + + 2004/04/28: Javier Rico * mifae/Changelog - Added. Next changes in directory mifae will be reported Index: trunk/MagicSoft/Mars/mtemp/mifae/Makefile =================================================================== --- trunk/MagicSoft/Mars/mtemp/mifae/Makefile (revision 3873) +++ trunk/MagicSoft/Mars/mtemp/mifae/Makefile (revision 3874) @@ -30,4 +30,5 @@ INCLUDES = -I. \ -I../../mbase \ + -I../../mjobs \ -I../../mpedestal \ -I../../mbadpixels \ Index: trunk/MagicSoft/Mars/mtemp/mifae/makeHillas.cc =================================================================== --- trunk/MagicSoft/Mars/mtemp/mifae/makeHillas.cc (revision 3873) +++ trunk/MagicSoft/Mars/mtemp/mifae/makeHillas.cc (revision 3874) @@ -19,5 +19,5 @@ #include "MHCalibrationChargeCam.h" #include "MHCalibrationRelTimeCam.h" -#include "MExtractSignal.h" +#include "MExtractFixedWindow.h" #include "MCalibrationChargeCalc.h" #include "MFCosmics.h" @@ -39,7 +39,11 @@ #include "MProgressBar.h" #include "MArgs.h" +#include "MRunIter.h" +#include "MJPedestal.h" +#include "MJCalibration.h" #include #include +#include using namespace std; @@ -56,7 +60,8 @@ // declaration of variables read from datacards TString outname; -TString pedname; -TString calname; -TString* datafile; +TString idirname; +MRunIter caliter; +MRunIter pediter; +MRunIter datiter; ULong_t nmaxevents=999999999; Short_t calflag=1; @@ -66,4 +71,7 @@ const TString defaultcard="input.datacard"; + +/*************************************************************/ +// makeHillas usage output static void Usage() { @@ -75,4 +83,6 @@ } +/*************************************************************/ +// main program int main(int argc, char **argv) { @@ -85,8 +95,10 @@ } + // get name of input datacard file TString datacard = arg.GetArgumentStr(0); if(!datacard.Length()) datacard = defaultcard; + // read the datacards if(!readDatacards(datacard)) { @@ -94,50 +106,45 @@ return -1; } + + // make the hillas file makeHillas(); } +/*************************************************************/ +/*************************************************************/ /*************************************************************/ void makeHillas() { -#if 0 - MStatusDisplay *display = new MStatusDisplay; - display->SetUpdateTime(3000); - display->Resize(850,700); - - gStyle->SetOptStat(1111); - gStyle->SetOptFit(); -#endif + + // set the signal extractor and calibration mode + MExtractFixedWindow extractor; + extractor.SetRange(hifirst,hilast,lofirst,lolast); + + MCalibrate::CalibrationMode_t calMode=MCalibrate::kDefault; + if(calflag==0) + calMode=MCalibrate::kNone; + if(calflag==-1) + calMode=MCalibrate::kDummy; + MCalibrate calibrate(calMode); + + + // general containers + MBadPixelsCam badcam; + MGeomCamMagic geomcam; + MGeomApply geomapl; + + // If you want to exclude pixels from the beginning, read + // an ascii-file with the corr. pixel numbers (see MBadPixelsCam) + // + // badcam.AsciiRead("badpixels.dat"); + /************************************/ /* FIRST LOOP: PEDESTAL COMPUTATION */ /************************************/ - - MParList plist1; - MTaskList tlist1; - plist1.AddToList(&tlist1); - - // containers - MPedestalCam pedcam; - MBadPixelsCam badcam; - // badcam.AsciiRead("badpixels.dat"); - - plist1.AddToList(&pedcam); - plist1.AddToList(&badcam); - - //tasks - MReadMarsFile read("Events", pedname); - read.DisableAutoScheme(); - - MGeomApply geomapl; - MPedCalcPedRun pedcalc; - - tlist1.AddToList(&read); - tlist1.AddToList(&geomapl); - tlist1.AddToList(&pedcalc); - - // Create and setup the eventloop - MEvtLoop pedloop; - pedloop.SetParList(&plist1); - //pedloop.SetDisplay(display); + + MJPedestal pedloop; + pedloop.SetInput(&pediter); + pedloop.SetBadPixels(badcam); cout << "*************************" << endl; @@ -145,332 +152,35 @@ cout << "*************************" << endl; - if (!pedloop.Eventloop()) + if (!pedloop.Process()) return; - - tlist1.PrintStatistics(); /*****************************/ /* SECOND LOOP: CALIBRATION */ /*****************************/ - - MParList plist2; - MTaskList tlist2; - plist2.AddToList(&tlist2); - plist2.AddToList(&pedcam); - plist2.AddToList(&badcam); - - // new containers - MGeomCamMagic geomcam; - MExtractedSignalCam sigcam; - MCalibrationChargeCam calcam; - MHCalibrationChargeCam histcharge; - MHCalibrationRelTimeCam histtime; - - plist2.AddToList(&geomcam); - plist2.AddToList(&sigcam); - plist2.AddToList(&calcam); - plist2.AddToList(&histcharge); - - //tasks - MReadMarsFile read2("Events", calname); - read2.DisableAutoScheme(); - - MExtractSignal sigcalc; - sigcalc.SetRange(hifirst,hilast,lofirst,lolast); - MCalibrationChargeCalc calcalc; - MFCosmics cosmics; - MContinue cont(&cosmics); - - - MFillH fillcam ("MHCalibrationChargeCam" , "MExtractedSignalCam"); - - tlist2.AddToList(&read2); - tlist2.AddToList(&geomapl); - tlist2.AddToList(&sigcalc); - tlist2.AddToList(&cont); - tlist2.AddToList(&fillcam); - tlist2.AddToList(&calcalc); - - // Create and setup the eventloop - MEvtLoop calloop; - calloop.SetParList(&plist2); - //calloop.SetDisplay(display); - + MCalibrationQECam qecam; + MJCalibration calloop; + calloop.SetInput(&caliter); + calloop.SetExtractor(&extractor); + // + // Set the corr. cams: + // + calloop.SetQECam(qecam); + calloop.SetBadPixels(pedloop.GetBadPixels()); + + // + // Apply rel. time calibration: + // calloop.SetRelTimeCalibration(); + + // Use as arrival time extractor MArrivalTimeCalc2: + // calloop.SetArrivalTimeLevel(2); + cout << "***************************" << endl; cout << "** COMPUTING CALIBRATION **" << endl; cout << "***************************" << endl; - if (!calloop.Eventloop()) + if (!calloop.Process(pedloop.GetPedestalCam())) return; - tlist2.PrintStatistics(); - - MLog gauglog; - gauglog.SetOutputFile(Form("%s%s",calcam.GetName(),".txt"),1); - calcam.SetLogStream(&gauglog); - calcam.Print(); - calcam.SetLogStream(&gLog); -#if 0 - // Create histograms to display - MHCamera disp1 (geomcam, "Cal;Charge", "Fitted Mean Charges"); - MHCamera disp2 (geomcam, "Cal;SigmaCharge", "Sigma of Fitted Charges"); - MHCamera disp3 (geomcam, "Cal;FitProb", "Probability of Fit"); - MHCamera disp4 (geomcam, "Cal;RSigma", "Reduced Sigmas"); - MHCamera disp5 (geomcam, "Cal;RSigma/Charge", "Reduced Sigma per Charge"); - MHCamera disp6 (geomcam, "Cal;FFactorPh", "Nr. of Photo-electrons (F-Factor Method)"); - MHCamera disp7 (geomcam, "Cal;FFactorConv", "Conversion Factor to photons (F-Factor Method)"); - MHCamera disp8 (geomcam, "Cal;FFactorFFactor", "Total F-Factor (F-Factor Method)"); - MHCamera disp9 (geomcam, "Cal;BlindPixPh", "Photon flux inside plexiglass (Blind Pixel Method)"); - MHCamera disp10 (geomcam, "Cal;BlindPixConv", "Conversion Factor to photons (Blind Pixel Method)"); - MHCamera disp11 (geomcam, "Cal;BlindPixFFactor","Total F-Factor (Blind Pixel Method)"); - MHCamera disp12 (geomcam, "Cal;PINDiodePh", "Photon flux outside plexiglass (PIN Diode Method)"); - MHCamera disp13 (geomcam, "Cal;PINDiodeConv", "Conversion Factor tp photons (PIN Diode Method)"); - MHCamera disp14 (geomcam, "Cal;PINDiodeFFactor","Total F-Factor (PIN Diode Method)"); - MHCamera disp15 (geomcam, "Cal;Excluded", "Pixels previously excluded"); - MHCamera disp16 (geomcam, "Cal;NotFitted", "Pixels that could not be fitted"); - MHCamera disp17 (geomcam, "Cal;NotFitValid", "Pixels with not valid fit results"); - MHCamera disp18 (geomcam, "Cal;HiGainOscillating", "Oscillating Pixels HI Gain"); - MHCamera disp19 (geomcam, "Cal;LoGainOscillating", "Oscillating Pixels LO Gain"); - MHCamera disp20 (geomcam, "Cal;HiGainPickup", "Number Pickup events Hi Gain"); - MHCamera disp21 (geomcam, "Cal;LoGainPickup", "Number Pickup events Lo Gain"); - MHCamera disp22 (geomcam, "Cal;Saturation", "Pixels with saturated Hi Gain"); - MHCamera disp23 (geomcam, "Cal;FFactorValid", "Pixels with valid F-Factor calibration"); - MHCamera disp24 (geomcam, "Cal;BlindPixelValid", "Pixels with valid BlindPixel calibration"); - MHCamera disp25 (geomcam, "Cal;PINdiodeFFactorValid", "Pixels with valid PINDiode calibration"); - - MHCamera disp26 (geomcam, "Cal;Ped", "Pedestals"); - MHCamera disp27 (geomcam, "Cal;PedRms", "Pedestal RMS"); - - MHCamera disp28 (geomcam, "time;Time", "Rel. Arrival Times"); - MHCamera disp29 (geomcam, "time;SigmaTime", "Sigma of Rel. Arrival Times"); - MHCamera disp30 (geomcam, "time;TimeProb", "Probability of Time Fit"); - MHCamera disp31 (geomcam, "time;NotFitValid", "Pixels with not valid fit results"); - MHCamera disp32 (geomcam, "time;Oscillating", "Oscillating Pixels"); - - MHCamera disp33 (geomcam, "Cal;AbsTimeMean", "Abs. Arrival Times"); - MHCamera disp34 (geomcam, "Cal;AbsTimeRms", "RMS of Arrival Times"); - - // Fitted charge means and sigmas - disp1.SetCamContent(calcam, 0); - disp1.SetCamError( calcam, 1); - disp2.SetCamContent(calcam, 2); - disp2.SetCamError( calcam, 3); - - // Fit probabilities - disp3.SetCamContent(calcam, 4); - - // Reduced Sigmas and reduced sigmas per charge - disp4.SetCamContent(calcam, 5); - disp4.SetCamError( calcam, 6); - disp5.SetCamContent(calcam, 7); - disp5.SetCamError( calcam, 8); - - // F-Factor Method - disp6.SetCamContent(calcam, 9); - disp6.SetCamError( calcam, 10); - disp7.SetCamContent(calcam, 11); - disp7.SetCamError( calcam, 12); - disp8.SetCamContent(calcam, 13); - disp8.SetCamError( calcam, 14); - - // Blind Pixel Method - disp9.SetCamContent(calcam, 15); - disp9.SetCamError( calcam, 16); - disp10.SetCamContent(calcam,17); - disp10.SetCamError( calcam,18); - disp11.SetCamContent(calcam,19); - disp11.SetCamError( calcam,20); - - // PIN Diode Method - disp12.SetCamContent(calcam,21); - disp12.SetCamError( calcam,22); - disp13.SetCamContent(calcam,23); - disp13.SetCamError( calcam,24); - disp14.SetCamContent(calcam,25); - disp14.SetCamError( calcam,26); - - // Pixels with defects - disp15.SetCamContent(calcam,27); - disp16.SetCamContent(calcam,28); - disp17.SetCamContent(badcam,9); - disp18.SetCamContent(badcam,16); - disp19.SetCamContent(badcam,15); - disp20.SetCamContent(calcam,29); - disp21.SetCamContent(calcam,30); - - // Lo Gain calibration - disp22.SetCamContent(calcam,31); - - // Valid flags - disp23.SetCamContent(calcam,32); - disp24.SetCamContent(calcam,33); - disp25.SetCamContent(calcam,34); - - // Pedestals - disp26.SetCamContent(calcam,35); - disp26.SetCamError( calcam,36); - disp27.SetCamContent(calcam,37); - disp27.SetCamError( calcam,38); - - // Relative Times - disp28.SetCamContent(histtime,0); - disp28.SetCamError( histtime,1); - disp29.SetCamContent(histtime,2); - disp29.SetCamError( histtime,3); - disp30.SetCamContent(histtime,4); - disp31.SetCamContent(histtime,5); - disp32.SetCamContent(histtime,6); - - // Absolute Times - disp33.SetCamContent(calcam,39); - disp33.SetCamError( calcam,40); - disp34.SetCamContent(calcam,41); - - disp1.SetYTitle("Charge [FADC units]"); - disp2.SetYTitle("\\sigma_{Charge} [FADC units]"); - disp3.SetYTitle("P_{Charge} [1]"); - - disp4.SetYTitle("\\sqrt{\\sigma^{2}_{Charge} - RMS^{2}_{Ped}} [FADC Counts]"); - disp5.SetYTitle("Reduced Sigma / Mean Charge [1]"); - - disp6.SetYTitle("Nr. Photo-electrons [1]"); - disp7.SetYTitle("Conversion Factor [Ph/FADC Count]"); - disp8.SetYTitle("\\sqrt{N_{Ph}}*\\sigma_{Charge}/\\mu_{Charge} [1] "); - - disp9.SetYTitle("Photon flux [ph/mm^2]"); - disp10.SetYTitle("Conversion Factor [Phot/FADC Count]"); - disp11.SetYTitle("\\sqrt{N_{Ph}}*\\sigma_{Charge}/\\mu_{Charge} [1]"); - - disp12.SetYTitle("Photon flux [ph/mm^2]"); - disp13.SetYTitle("Conversion Factor [Phot/FADC Count]"); - disp14.SetYTitle("\\sqrt{N_{Ph}}*\\sigma_{Charge}/\\mu_{Charge} [1]"); - - disp15.SetYTitle("[1]"); - disp16.SetYTitle("[1]"); - disp17.SetYTitle("[1]"); - disp18.SetYTitle("[1]"); - disp19.SetYTitle("[1]"); - disp20.SetYTitle("[1]"); - disp21.SetYTitle("[1]"); - disp22.SetYTitle("[1]"); - disp23.SetYTitle("[1]"); - disp24.SetYTitle("[1]"); - disp25.SetYTitle("[1]"); - - disp26.SetYTitle("Ped [FADC Counts ]"); - disp27.SetYTitle("RMS_{Ped} [FADC Counts ]"); - - disp28.SetYTitle("Time Offset [ns]"); - disp29.SetYTitle("Timing resolution [ns]"); - disp30.SetYTitle("P_{Time} [1]"); - - disp31.SetYTitle("[1]"); - disp32.SetYTitle("[1]"); - - disp33.SetYTitle("Mean Abs. Time [FADC slice]"); - disp34.SetYTitle("RMS Abs. Time [FADC slices]"); - - gStyle->SetOptStat(1111); - gStyle->SetOptFit(); - - // Charges - TCanvas &c1 = display->AddTab("Fit.Charge"); - c1.Divide(2, 3); - - CamDraw(c1, disp1,calcam,1, 2 , 2); - CamDraw(c1, disp2,calcam,2, 2 , 2); - - // Fit Probability - TCanvas &c2 = display->AddTab("Fit.Prob"); - c2.Divide(1,3); - - CamDraw(c2, disp3,calcam,1, 1 , 4); - - // Reduced Sigmas - TCanvas &c3 = display->AddTab("Red.Sigma"); - c3.Divide(2,3); - - CamDraw(c3, disp4,calcam,1, 2 , 2); - CamDraw(c3, disp5,calcam,2, 2 , 2); - - // F-Factor Method - TCanvas &c4 = display->AddTab("F-Factor"); - c4.Divide(3,3); - - CamDraw(c4, disp6,calcam,1, 3 , 2); - CamDraw(c4, disp7,calcam,2, 3 , 2); - CamDraw(c4, disp8,calcam,3, 3 , 2); - - // Blind Pixel Method - TCanvas &c5 = display->AddTab("BlindPix"); - c5.Divide(3, 3); - - CamDraw(c5, disp9,calcam,1, 3 , 9); - CamDraw(c5, disp10,calcam,2, 3 , 2); - CamDraw(c5, disp11,calcam,3, 3 , 2); - - // PIN Diode Method - TCanvas &c6 = display->AddTab("PINDiode"); - c6.Divide(3,3); - - CamDraw(c6, disp12,calcam,1, 3 , 9); - CamDraw(c6, disp13,calcam,2, 3 , 2); - CamDraw(c6, disp14,calcam,3, 3 , 2); - - // Defects - TCanvas &c7 = display->AddTab("Defects"); - c7.Divide(4,2); - - CamDraw(c7, disp15,calcam,1,4, 0); - CamDraw(c7, disp16,calcam,2,4, 0); - CamDraw(c7, disp20,calcam,3,4, 0); - CamDraw(c7, disp21,calcam,4,4, 0); - - // BadCam - TCanvas &c8 = display->AddTab("Defects"); - c8.Divide(3,2); - - CamDraw(c8, disp17,badcam,1,3, 0); - CamDraw(c8, disp18,badcam,2,3, 0); - CamDraw(c8, disp19,badcam,3,3, 0); - - // Valid flags - TCanvas &c9 = display->AddTab("Validity"); - c9.Divide(4,2); - - CamDraw(c9, disp22,calcam,1,4,0); - CamDraw(c9, disp23,calcam,2,4,0); - CamDraw(c9, disp24,calcam,3,4,0); - CamDraw(c9, disp25,calcam,4,4,0); - - // Pedestals - TCanvas &c10 = display->AddTab("Pedestals"); - c10.Divide(2,3); - - CamDraw(c10,disp26,calcam,1,2,1); - CamDraw(c10,disp27,calcam,2,2,2); - - // Rel. Times - TCanvas &c11 = display->AddTab("Fitted Rel. Times"); - c11.Divide(3,3); - - CamDraw(c11,disp28,calcam,1,3,2); - CamDraw(c11,disp29,calcam,2,3,2); - CamDraw(c11,disp30,calcam,3,3,4); - - // Time Defects - TCanvas &c12 = display->AddTab("Time Def."); - c12.Divide(2,2); - - CamDraw(c12, disp31,calcam,1,2, 0); - CamDraw(c12, disp32,calcam,2,2, 0); - - // Abs. Times - TCanvas &c13 = display->AddTab("Abs. Times"); - c13.Divide(2,3); - - CamDraw(c13,disp33,calcam,1,2,2); - CamDraw(c13,disp34,calcam,2,2,2); -#endif + /************************************************************************/ /* THIRD LOOP: PEDESTAL CALIBRATION INTO PHOTONS */ @@ -485,30 +195,23 @@ plist3.AddToList(&geomcam); - plist3.AddToList(&pedcam); - plist3.AddToList(&calcam); - plist3.AddToList(&badcam); - plist3.AddToList(&sigcam); + plist3.AddToList(&pedloop.GetPedestalCam()); + plist3.AddToList(&calloop.GetCalibrationCam()); + plist3.AddToList(&calloop.GetQECam()); + plist3.AddToList(&calloop.GetRelTimeCam()); + plist3.AddToList(&calloop.GetBadPixels()); plist3.AddToList(&nphot); plist3.AddToList(&nphotrms); - - MCalibrate::CalibrationMode_t calMode=MCalibrate::kDefault; - if(calflag==0) - calMode=MCalibrate::kNone; - - //tasks - MReadMarsFile read3("Events", pedname); + MReadMarsFile read3("Events"); read3.DisableAutoScheme(); - - MExtractSignal extsig; - extsig.SetRange(hifirst,hilast,lofirst,lolast); - MCalibrate photcalc(calMode); + static_cast(read3).AddFiles(pediter); + MPedPhotCalc photrmscalc; tlist3.AddToList(&read3); tlist3.AddToList(&geomapl); - tlist3.AddToList(&extsig); - tlist3.AddToList(&photcalc); + tlist3.AddToList(&extractor); + tlist3.AddToList(&calibrate); tlist3.AddToList(&photrmscalc); @@ -531,14 +234,16 @@ // containers MHillas hillas; - MSrcPosCam source; + MSrcPosCam srcposcam; MRawRunHeader runhead; - + plist4.AddToList(&geomcam); - plist4.AddToList(&pedcam); - plist4.AddToList(&calcam); - plist4.AddToList(&badcam); + plist4.AddToList(&pedloop.GetPedestalCam()); + plist4.AddToList(&calloop.GetCalibrationCam()); + plist4.AddToList(&calloop.GetQECam()); + plist4.AddToList(&calloop.GetRelTimeCam()); + plist4.AddToList(&calloop.GetBadPixels()); plist4.AddToList(&nphot); plist4.AddToList(&nphotrms); - plist4.AddToList(&source); + plist4.AddToList(&srcposcam); plist4.AddToList(&hillas); plist4.AddToList(&runhead); @@ -546,19 +251,6 @@ //tasks MReadMarsFile read4("Events"); - - for(Int_t i=0;i(read4).AddFiles(datiter); MImgCleanStd clean(lcore,ltail); @@ -568,34 +260,19 @@ MWriteRootFile write(outname,"RECREATE"); - // write.AddContainer("MGeomCam" , "RunHeaders"); - // write.AddContainer("MRawRunHeader" , "RunHeaders"); - // write.AddContainer("MSrcPosCam" , "RunHeaders"); - // write.AddContainer("MCalibrationChargeCam" , "RunHeaders"); - // write.AddContainer("MPedPhotCam" , "RunHeaders"); // Attention, was in Events - Tree!! - // write.AddContainer("MPedestalCam" , "RunHeaders"); - // write.AddContainer("MHCalibrationRelTimeCam","RunHeaders"); - - // write.AddContainer("MCerPhotEvt" , "Events"); - // write.AddContainer("MRawEvtHeader" , "Events"); - // write.AddContainer("MBadPixelsCam" , "Events"); - // write.AddContainer("MPedPhotCam" , "Events"); - - write.AddContainer("MHillas" , "Parameters"); - write.AddContainer("MHillasSrc" , "Parameters"); - write.AddContainer("MHillasExt" , "Parameters"); - write.AddContainer("MNewImagePar" , "Parameters"); - write.AddContainer("MRawEvtHeader" , "Parameters"); - write.AddContainer("MRawRunHeader" , "Parameters"); + write.AddContainer("MHillas" , "Parameters"); + write.AddContainer("MHillasSrc" , "Parameters"); + write.AddContainer("MHillasExt" , "Parameters"); + write.AddContainer("MNewImagePar" , "Parameters"); + write.AddContainer("MRawEvtHeader" , "Parameters"); + write.AddContainer("MRawRunHeader" , "Parameters"); write.AddContainer("MConcentration" , "Parameters"); + write.AddContainer("MSrcPosCam", "Parameters"); tlist4.AddToList(&read4); tlist4.AddToList(&geomapl); - tlist4.AddToList(&extsig); - tlist4.AddToList(&photcalc); - //tlist4.AddToList(&blind); + tlist4.AddToList(&extractor); + tlist4.AddToList(&calibrate); tlist4.AddToList(&clean); - //tlist4.AddToList(&blind2); tlist4.AddToList(&hcalc); - // tlist4.AddToList(&srcposcalc); tlist4.AddToList(&csrc1); tlist4.AddToList(&write); @@ -604,6 +281,4 @@ MEvtLoop datloop; datloop.SetParList(&plist4); - // MProgressBar bar; - // datloop.SetProgressBar(&bar); cout << "*************************************************************" << endl; @@ -615,152 +290,7 @@ tlist4.PrintStatistics(); - } -#if 0 -void CamDraw(TCanvas &c, MHCamera &cam, MCamEvent &evt, Int_t i, Int_t j, Int_t fit) -{ - - c.cd(i); - gPad->SetBorderMode(0); - MHCamera *obj1=(MHCamera*)cam.DrawCopy("hist"); - // obj1->AddNotify(evt); - - c.cd(i+j); - gPad->SetBorderMode(0); - obj1->Draw(); - ((MHCamera*)obj1)->SetPrettyPalette(); - - if (fit != 0) - { - c.cd(i+2*j); - gPad->SetBorderMode(0); - TH1D *obj2 = (TH1D*)obj1->Projection(obj1.GetName()); - -// obj2->Sumw2(); - obj2->Draw(); - obj2->SetBit(kCanDelete); - - const Double_t min = obj2->GetBinCenter(obj2->GetXaxis()->GetFirst()); - const Double_t max = obj2->GetBinCenter(obj2->GetXaxis()->GetLast()); - const Double_t integ = obj2->Integral("width")/2.5066283; - const Double_t mean = obj2->GetMean(); - const Double_t rms = obj2->GetRMS(); - const Double_t width = max-min; - - if (rms == 0. || width == 0. ) - return; - - switch (fit) - { - case 1: - TF1 *sgaus = new TF1("sgaus","gaus(0)",min,max); - sgaus->SetBit(kCanDelete); - sgaus->SetParNames("Area","#mu","#sigma"); - sgaus->SetParameters(integ/rms,mean,rms); - sgaus->SetParLimits(0,0.,integ); - sgaus->SetParLimits(1,min,max); - sgaus->SetParLimits(2,0,width/1.5); - obj2->Fit("sgaus","QLR"); - obj2->GetFunction("sgaus")->SetLineColor(kYellow); - break; - - case 2: - TString dgausform = "([0]-[3])/[2]*exp(-0.5*(x-[1])*(x-[1])/[2]/[2])"; - dgausform += "+[3]/[5]*exp(-0.5*(x-[4])*(x-[4])/[5]/[5])"; - TF1 *dgaus = new TF1("dgaus",dgausform.Data(),min,max); - dgaus->SetBit(kCanDelete); - dgaus->SetParNames("A_{tot}","#mu_{1}","#sigma_{1}","A_{2}","#mu_{2}","#sigma_{2}"); - dgaus->SetParameters(integ,(min+mean)/2.,width/4., - integ/width/2.,(max+mean)/2.,width/4.); - // The left-sided Gauss - dgaus->SetParLimits(0,integ-1.5,integ+1.5); - dgaus->SetParLimits(1,min+(width/10.),mean); - dgaus->SetParLimits(2,0,width/2.); - // The right-sided Gauss - dgaus->SetParLimits(3,0,integ); - dgaus->SetParLimits(4,mean,max-(width/10.)); - dgaus->SetParLimits(5,0,width/2.); - obj2->Fit("dgaus","QLRM"); - obj2->GetFunction("dgaus")->SetLineColor(kYellow); - break; - - case 3: - TString tgausform = "([0]-[3]-[6])/[2]*exp(-0.5*(x-[1])*(x-[1])/[2]/[2])"; - tgausform += "+[3]/[5]*exp(-0.5*(x-[4])*(x-[4])/[5]/[5])"; - tgausform += "+[6]/[8]*exp(-0.5*(x-[7])*(x-[7])/[8]/[8])"; - TF1 *tgaus = new TF1("tgaus",tgausform.Data(),min,max); - tgaus->SetBit(kCanDelete); - tgaus->SetParNames("A_{tot}","#mu_{1}","#sigma_{1}", - "A_{2}","#mu_{2}","#sigma_{2}", - "A_{3}","#mu_{3}","#sigma_{3}"); - tgaus->SetParameters(integ,(min+mean)/2,width/4., - integ/width/3.,(max+mean)/2.,width/4., - integ/width/3.,mean,width/2.); - // The left-sided Gauss - tgaus->SetParLimits(0,integ-1.5,integ+1.5); - tgaus->SetParLimits(1,min+(width/10.),mean); - tgaus->SetParLimits(2,width/15.,width/2.); - // The right-sided Gauss - tgaus->SetParLimits(3,0.,integ); - tgaus->SetParLimits(4,mean,max-(width/10.)); - tgaus->SetParLimits(5,width/15.,width/2.); - // The Gauss describing the outliers - tgaus->SetParLimits(6,0.,integ); - tgaus->SetParLimits(7,min,max); - tgaus->SetParLimits(8,width/4.,width/1.5); - obj2->Fit("tgaus","QLRM"); - obj2->GetFunction("tgaus")->SetLineColor(kYellow); - break; - case 4: - obj2->Fit("pol0","Q"); - obj2->GetFunction("pol0")->SetLineColor(kYellow); - break; - case 9: - break; - default: - obj2->Fit("gaus","Q"); - obj2->GetFunction("gaus")->SetLineColor(kYellow); - break; - } - - TArrayI s0(3); - s0[0] = 6; - s0[1] = 1; - s0[2] = 2; - - TArrayI s1(3); - s1[0] = 3; - s1[1] = 4; - s1[2] = 5; - - TArrayI inner(1); - inner[0] = 0; - - TArrayI outer(1); - outer[0] = 1; - - // Just to get the right (maximum) binning - TH1D *half[4]; - half[0] = obj1->ProjectionS(s0, inner, "Sector 6-1-2 Inner"); - half[1] = obj1->ProjectionS(s1, inner, "Sector 3-4-5 Inner"); - half[2] = obj1->ProjectionS(s0, outer, "Sector 6-1-2 Outer"); - half[3] = obj1->ProjectionS(s1, outer, "Sector 3-4-5 Outer"); - - for (int i=0; i<4; i++) - { - half[i]->SetLineColor(kRed+i); - half[i]->SetDirectory(0); - half[i]->SetBit(kCanDelete); - half[i]->Draw("same"); - } - - gPad->Modified(); - gPad->Update(); - - } -} -#endif - +/******************************************************************************/ Bool_t readDatacards(TString& filename) { @@ -773,5 +303,4 @@ TString word; - Int_t ifile = 0; while(ifun >> word) @@ -788,49 +317,38 @@ ifun >> nmaxevents; - // pedestal file name - if(strcmp(word.Data(),"PFILE")==0) - { - if(pedname.Length()) - cout << "readDataCards Warning: overriding pedestal file name" << endl; - ifun >> pedname; - } - - // calibration file name - if(strcmp(word.Data(),"CFILE")==0) - { - if(calname.Length()) - cout << "readDataCards Warning: overriding calibration file name" << endl; - ifun >> calname; - } - - // number of data files - if(strcmp(word.Data(),"NFILES")==0) - { - if(nfiles) - { - Int_t dummy; - cout << "readDataCards Warning: trying to set a new value of number of files. Skiping..." << endl; - ifun >> dummy; - } - else - { - ifun >> nfiles; - datafile = new TString[nfiles]; - } - } - - - // data file name - if(strcmp(word.Data(),"DFILE")==0) - { - TString newfile; - ifun >> newfile; - if(ifile> idirname; + } + + // pedestal runs + if(strcmp(word.Data(),"PRUNS")==0) + { + if(pediter.GetNumRuns()) + cout << "readDataCards Warning: adding pedestal runs to the existing list" << endl; + ifun >> word; + pediter.AddRuns(word.Data(),idirname.Data()); + } + + // calibration runs + if(strcmp(word.Data(),"CRUNS")==0) + { + if(caliter.GetNumRuns()) + cout << "readDataCards Warning: adding calibration runs to the existing list" << endl; + ifun >> word; + caliter.AddRuns(word.Data(),idirname.Data()); + } + + // data runs + if(strcmp(word.Data(),"DRUNS")==0) + { + if(datiter.GetNumRuns()) + cout << "readDataCards Warning: adding data runs to the existing list" << endl; + ifun >> word; + datiter.AddRuns(word.Data(),idirname.Data()); } @@ -855,4 +373,8 @@ } + pediter.Reset(); + caliter.Reset(); + datiter.Reset(); + TString pfile; // Dump read values @@ -860,10 +382,13 @@ cout << "* Datacards read from file " << filename << endl; cout << "************************************************" << endl; - cout << "Pedestal file: " << pedname << endl; - cout << "Calibration file: " << calname << endl; - cout << "Number of data files: " << nfiles << endl; - cout << "Data files: " << endl; - for(int i=0;i