Changeset 4609
- Timestamp:
- 08/13/04 14:59:17 (20 years ago)
- Location:
- trunk/MagicSoft/Mars
- Files:
-
- 22 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/MagicSoft/Mars/macros/readraw.C
r948 r4609 16 16 ! 17 17 ! 18 ! Author(s): Thomas Bretz 12/2000 (tbretz@uni-sw.gwdg.de)18 ! Author(s): Thomas Bretz 12/2000 <mailto:tbretz@astro.uni-wuerzburg.de> 19 19 ! 20 ! Copyright: MAGIC Software Development, 2000-200 120 ! Copyright: MAGIC Software Development, 2000-2004 21 21 ! 22 22 ! … … 24 24 25 25 26 void readraw( )26 void readraw(const char *fname="/data/MAGIC/Period016/mcdata/spot_1cm/standard/gamma/Gamma_zbin9_90_7_1740to1749_w0.root") 27 27 { 28 28 // 29 29 // open the file 30 30 // 31 TFile input( "test.root", "READ");31 TFile input(fname, "READ"); 32 32 33 33 // … … 46 46 // 47 47 MRawRunHeader *runheader = new MRawRunHeader(); 48 runtree->GetBranch("MRawRunHeader ")->SetAddress(&runheader);48 runtree->GetBranch("MRawRunHeader.")->SetAddress(&runheader); 49 49 runtree->GetEvent(0); 50 50 runheader->Print(); … … 53 53 // open the Data Tree 54 54 // 55 TTree *evttree = (TTree*) input.Get(" Data");55 TTree *evttree = (TTree*) input.Get("Events"); 56 56 57 57 // 58 58 // create the instances of the data to read in 59 59 // 60 MRawEvtHeader *evtheader = new MRawEvtHeader();61 MTime *evttime = new MTime();62 MRawEvtData *evtdata = new MRawEvtData();63 MRawCrateArray *evtcrate = new MRawCrateArray();60 MRawEvtHeader *evtheader = 0; 61 MTime *evttime = 0; 62 MRawEvtData *evtdata = 0; 63 MRawCrateArray *evtcrate = 0; 64 64 65 65 // 66 66 // enable the corresponding branches 67 67 // 68 evttree->GetBranch("MRawEvtHeader")->SetAddress(&evtheader); 69 evttree->GetBranch("MTime")->SetAddress(&evttime); 70 evttree->GetBranch("MRawEvtData")->SetAddress(&evtdata); 71 evttree->GetBranch("MRawCrateArray")->SetAddress(&evtcrate); 68 evttree->GetBranch("MRawEvtHeader.")->SetAddress(&evtheader); 69 evttree->GetBranch("MRawEvtData.")->SetAddress(&evtdata); 70 71 // Use this for real data only 72 //evttree->GetBranch("MTime.")->SetAddress(&evttime); 73 //evttree->GetBranch("MRawCrateArray.")->SetAddress(&evtcrate); 72 74 73 75 // … … 83 85 84 86 evtheader->Print(); 85 evttime->Print();86 evtcrate->Print();87 87 evtdata->Print(); 88 89 // Use this for real data only 90 //evttime->Print(); 91 //evtcrate->Print(); 88 92 } 89 93 } -
trunk/MagicSoft/Mars/mbase/MTaskEnv.cc
r4601 r4609 172 172 if (!fTask) 173 173 { 174 *fLog << err << GetDescriptor() << " - ERROR: No ttask setup." << endl;174 *fLog << err << GetDescriptor() << " - ERROR: No task setup." << endl; 175 175 return kFALSE; 176 176 } -
trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCalc.cc
r4603 r4609 1524 1524 } 1525 1525 1526 *fLog << inf << GetDescriptor() << ": Mean overall F-Factor (photons to FADC counts) " 1527 << "with area index: " << i << ": " 1528 << Form("%4.2f+-%4.2f",mean,sigma) << endl; 1526 *fLog << inf << endl; 1527 *fLog << inf << GetDescriptor() << ": Mean F-Factor " 1528 << "with area index #" << i << ": " 1529 << Form("%4.2f+-%4.2f",mean,sigma) << "fadc/ph" << endl; 1529 1530 1530 1531 lowlim [i] = 1.1; … … 1550 1551 if ( ffactor < lowlim[aidx] || ffactor > upplim[aidx] ) 1551 1552 { 1552 *fLog << warn << GetDescriptor() << ": Overall F-Factor: " 1553 << Form("%5.2f out of %3.1f sigma limit: ",ffactor,fFFactorErrLimit) 1554 << Form("[%5.2f,%5.2f] pixel %4i",lowlim[aidx],upplim[aidx],i) << endl; 1553 *fLog << warn << GetDescriptor() << ": Overall F-Factor " 1554 << Form("%5.2f",ffactor) << " out of range [" 1555 << Form("%5.2f,%5.2f",lowlim[aidx],upplim[aidx]) << "] pixel " << i << endl; 1556 1555 1557 bad.SetUncalibrated( MBadPixelsPix::kDeviatingFFactor ); 1556 1558 bad.SetUnsuitable ( MBadPixelsPix::kUnsuitableRun ); … … 1857 1859 return Form("%s/%s", (const char*)fOutputPath, (const char*)fOutputFile); 1858 1860 } 1861 1862 Int_t MCalibrationChargeCalc::ReadEnv(const TEnv &env, TString prefix, Bool_t print) 1863 { 1864 Bool_t rc = kFALSE; 1865 if (IsEnvDefined(env, prefix, "ChargeLimit", print)) 1866 { 1867 SetChargeLimit(GetEnvValue(env, prefix, "ChargeLimit", fChargeLimit)); 1868 rc = kTRUE; 1869 } 1870 1871 if (IsEnvDefined(env, prefix, "ChargeErrLimit", print)) 1872 { 1873 SetChargeErrLimit(GetEnvValue(env, prefix, "ChargeErrLimit", fChargeErrLimit)); 1874 rc = kTRUE; 1875 } 1876 1877 if (IsEnvDefined(env, prefix, "ChargeRelErrLimit", print)) 1878 { 1879 SetChargeRelErrLimit(GetEnvValue(env, prefix, "ChargeRelErrLimit", fChargeRelErrLimit)); 1880 rc = kTRUE; 1881 } 1882 if (IsEnvDefined(env, prefix, "Debug", print)) 1883 { 1884 SetDebug(GetEnvValue(env, prefix, "Debug", IsDebug())); 1885 rc = kTRUE; 1886 } 1887 if (IsEnvDefined(env, prefix, "FFactorErrLimit", print)) 1888 { 1889 SetFFactorErrLimit(GetEnvValue(env, prefix, "FFactorErrLimit", fFFactorErrLimit)); 1890 rc = kTRUE; 1891 } 1892 if (IsEnvDefined(env, prefix, "LambdaErrLimit", print)) 1893 { 1894 SetLambdaErrLimit(GetEnvValue(env, prefix, "LambdaErrLimit", fLambdaErrLimit)); 1895 rc = kTRUE; 1896 } 1897 if (IsEnvDefined(env, prefix, "LambdaCheckLimit", print)) 1898 { 1899 SetLambdaCheckLimit(GetEnvValue(env, prefix, "LambdaCheckLimit", fLambdaCheckLimit)); 1900 rc = kTRUE; 1901 } 1902 if (IsEnvDefined(env, prefix, "PheErrLimit", print)) 1903 { 1904 SetPheErrLimit(GetEnvValue(env, prefix, "PheErrLimit", fPheErrLimit)); 1905 rc = kTRUE; 1906 } 1907 // void SetPulserColor(const MCalibrationCam::PulserColor_t col) { fPulserColor = col; } 1908 1909 return rc; 1910 } -
trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCalc.h
r4603 r4609 127 127 Int_t Process (); 128 128 Int_t PostProcess(); 129 Int_t ReadEnv(const TEnv &env, TString prefix, Bool_t print); 129 130 130 131 public: -
trunk/MagicSoft/Mars/mcalib/MHCalibrationCam.cc
r4603 r4609 1023 1023 } 1024 1024 1025 Int_t MHCalibrationCam::ReadEnv(const TEnv &env, TString prefix, Bool_t print) 1026 { 1027 Bool_t rc = kFALSE; 1028 if (IsEnvDefined(env, prefix, "Debug", print)) 1029 { 1030 SetDebug(GetEnvValue(env, prefix, "Debug", IsDebug())); 1031 rc = kTRUE; 1032 } 1033 1034 return rc; 1035 } -
trunk/MagicSoft/Mars/mcalib/MHCalibrationCam.h
r4602 r4609 109 109 void InitHists(MHGausEvents &hist, MBadPixelsPix &bad, const Int_t i); 110 110 111 Int_t ReadEnv(const TEnv &env, TString prefix, Bool_t print); 112 111 113 public: 112 114 -
trunk/MagicSoft/Mars/mjobs/MJCalibration.cc
r4607 r4609 87 87 #include "MJCalibration.h" 88 88 89 #include <TEnv.h> 89 90 #include <TFile.h> 90 91 #include <TStyle.h> … … 96 97 97 98 #include "MRunIter.h" 99 #include "MSequence.h" 98 100 #include "MParList.h" 99 101 #include "MTaskList.h" … … 125 127 #include "MRawFileRead.h" 126 128 #include "MGeomApply.h" 129 #include "MTaskEnv.h" 127 130 #include "MBadPixelsMerge.h" 128 131 #include "MBadPixelsCam.h" … … 158 161 // 159 162 MJCalibration::MJCalibration(const char *name, const char *title) 160 : f Runs(0), fExtractor(NULL), fTimeExtractor(NULL),163 : fEnv(0), fRuns(0), fSequence(0), fExtractor(NULL), fTimeExtractor(NULL), 161 164 fColor(MCalibrationCam::kNONE), fDisplayType(kNormalDisplay), 162 165 fRelTimes(kFALSE), fDataCheck(kFALSE), fDebug(kFALSE) … … 171 174 } 172 175 176 MJCalibration::~MJCalibration() 177 { 178 if (fEnv) 179 delete fEnv; 180 } 173 181 174 182 // -------------------------------------------------------------------------- … … 224 232 TString title = fDisplay->GetTitle(); 225 233 title += "-- Calibration "; 226 title += f Runs->GetRunsAsString();234 title += fSequence ? Form("calib%06d", fSequence->GetSequence()) : fRuns->GetRunsAsString(); 227 235 title += " --"; 228 236 fDisplay->SetTitle(title); … … 234 242 235 243 // Create histograms to display 236 MHCamera disp1 (geomcam, Form("%s%s","Charge",(fRuns->GetRunsAsFileName()).Data()), 237 "Fitted Mean Charges"); 238 MHCamera disp2 (geomcam, Form("%s%s","SigmaCharge",(fRuns->GetRunsAsFileName()).Data()), 239 "Sigma of Fitted Charges"); 240 MHCamera disp3 (geomcam, Form("%s%s","RSigma",(fRuns->GetRunsAsFileName()).Data()), 241 "Reduced Sigmas"); 242 MHCamera disp4 (geomcam, Form("%s%s","RSigmaPerCharge",(fRuns->GetRunsAsFileName()).Data()), 243 "Reduced Sigma per Charge"); 244 MHCamera disp5 (geomcam, Form("%s%s","NumPhes",(fRuns->GetRunsAsFileName()).Data()), 245 "Nr. of Phe's (F-Factor Method)"); 246 MHCamera disp6 (geomcam, Form("%s%s","ConvFADC2Phes",(fRuns->GetRunsAsFileName()).Data()), 247 "Conversion Factor (F-Factor Method)"); 248 MHCamera disp7 (geomcam, Form("%s%s","TotalFFactor",(fRuns->GetRunsAsFileName()).Data()), 249 "Total F-Factor (F-Factor Method)"); 250 MHCamera disp8 (geomcam, Form("%s%s","CascadesQEFFactor",(fRuns->GetRunsAsFileName()).Data()), 251 "Cascades QE (F-Factor Method)"); 252 MHCamera disp9 (geomcam, Form("%s%s","CascadesQEBlindPix",(fRuns->GetRunsAsFileName()).Data()), 253 "Cascades QE (Blind Pixel Method)"); 254 MHCamera disp10(geomcam, Form("%s%s","CascadesQEPINDiode",(fRuns->GetRunsAsFileName()).Data()), 255 "Cascades QE (PIN Diode Method)"); 256 MHCamera disp11(geomcam, Form("%s%s","CascadesQECombined",(fRuns->GetRunsAsFileName()).Data()), 257 "Cascades QE (Combined Method)"); 258 MHCamera disp12(geomcam, Form("%s%s","FFactorValid",(fRuns->GetRunsAsFileName()).Data()), 259 "Pixels with valid F-Factor calibration"); 260 MHCamera disp13(geomcam, Form("%s%s","BlindPixelValid",(fRuns->GetRunsAsFileName()).Data()), 261 "Pixels with valid BlindPixel calibration"); 262 MHCamera disp14(geomcam, Form("%s%s","PINdiodeValid",(fRuns->GetRunsAsFileName()).Data()), 263 "Pixels with valid PINDiode calibration"); 264 MHCamera disp15(geomcam, Form("%s%s","CombinedValid",(fRuns->GetRunsAsFileName()).Data()), 265 "Pixels with valid Combined calibration"); 266 MHCamera disp16(geomcam, Form("%s%s","Saturation",(fRuns->GetRunsAsFileName()).Data()), 267 "Pixels with saturated Hi Gain"); 268 MHCamera disp17(geomcam, Form("%s%s","ConversionMeans",(fRuns->GetRunsAsFileName()).Data()), 269 "Conversion HiGain.vs.LoGain Means"); 270 MHCamera disp18(geomcam, Form("%s%s","ConversionSigmas",(fRuns->GetRunsAsFileName()).Data()), 271 "Conversion HiGain.vs.LoGain Sigmas"); 272 MHCamera disp19(geomcam, Form("%s%s","HiGainPickup",(fRuns->GetRunsAsFileName()).Data()), 273 "Number Pickup events Hi Gain"); 274 MHCamera disp20(geomcam, Form("%s%s","LoGainPickup",(fRuns->GetRunsAsFileName()).Data()), 275 "Number Pickup events Lo Gain"); 276 MHCamera disp21(geomcam, Form("%s%s","HiGainBlackout",(fRuns->GetRunsAsFileName()).Data()), 277 "Number Blackout events Hi Gain"); 278 MHCamera disp22(geomcam, Form("%s%s","LoGainBlackout",(fRuns->GetRunsAsFileName()).Data()), 279 "Number Blackout events Lo Gain"); 280 MHCamera disp23(geomcam, Form("%s%s","Excluded",(fRuns->GetRunsAsFileName()).Data()), 281 "Pixels previously excluded"); 282 MHCamera disp24(geomcam, Form("%s%s","UnSuitable",(fRuns->GetRunsAsFileName()).Data()), 283 "Pixels not suited for further analysis"); 284 MHCamera disp25(geomcam, Form("%s%s","UnReliable",(fRuns->GetRunsAsFileName()).Data()), 285 "Pixels not reliable for further analysis"); 286 MHCamera disp26(geomcam, Form("%s%s","HiGainOscillating",(fRuns->GetRunsAsFileName()).Data()), 287 "Oscillating Pixels High Gain"); 288 MHCamera disp27(geomcam, Form("%s%s","LoGainOscillating",(fRuns->GetRunsAsFileName()).Data()), 289 "Oscillating Pixels Low Gain"); 290 MHCamera disp28(geomcam, Form("%s%s","AbsTimeMean",(fRuns->GetRunsAsFileName()).Data()), 291 "Abs. Arrival Times"); 292 MHCamera disp29(geomcam, Form("%s%s","AbsTimeRms",(fRuns->GetRunsAsFileName()).Data()), 293 "RMS of Arrival Times"); 294 MHCamera disp30(geomcam, Form("%s%s","MeanTime",(fRuns->GetRunsAsFileName()).Data()), 295 "Mean Rel. Arrival Times"); 296 MHCamera disp31(geomcam, Form("%s%s","SigmaTime",(fRuns->GetRunsAsFileName()).Data()), 297 "Sigma Rel. Arrival Times"); 298 MHCamera disp32(geomcam, Form("%s%s","TimeProb",(fRuns->GetRunsAsFileName()).Data()), 299 "Probability of Time Fit"); 300 MHCamera disp33(geomcam, Form("%s%s","TimeNotFitValid",(fRuns->GetRunsAsFileName()).Data()), 301 "Pixels with not valid fit results"); 302 MHCamera disp34(geomcam, Form("%s%s","TimeOscillating",(fRuns->GetRunsAsFileName()).Data()), 303 "Oscillating Pixels"); 244 MHCamera disp1 (geomcam, "Charge", "Fitted Mean Charges"); 245 MHCamera disp2 (geomcam, "SigmaCharge", "Sigma of Fitted Charges"); 246 MHCamera disp3 (geomcam, "RSigma", "Reduced Sigmas"); 247 MHCamera disp4 (geomcam, "RSigmaPerCharge", "Reduced Sigma per Charge"); 248 MHCamera disp5 (geomcam, "NumPhes", "Nr. of Phe's (F-Factor Method)"); 249 MHCamera disp6 (geomcam, "ConvFADC2Phes", "Conversion Factor (F-Factor Method)"); 250 MHCamera disp7 (geomcam, "TotalFFactor", "Total F-Factor (F-Factor Method)"); 251 MHCamera disp8 (geomcam, "CascadesQEFFactor", "Cascades QE (F-Factor Method)"); 252 MHCamera disp9 (geomcam, "CascadesQEBlindPix","Cascades QE (Blind Pixel Method)"); 253 MHCamera disp10(geomcam, "CascadesQEPINDiode","Cascades QE (PIN Diode Method)"); 254 MHCamera disp11(geomcam, "CascadesQECombined","Cascades QE (Combined Method)"); 255 MHCamera disp12(geomcam, "FFactorValid", "Pixels with valid F-Factor calibration"); 256 MHCamera disp13(geomcam, "BlindPixelValid", "Pixels with valid BlindPixel calibration"); 257 MHCamera disp14(geomcam, "PINdiodeValid", "Pixels with valid PINDiode calibration"); 258 MHCamera disp15(geomcam, "CombinedValid", "Pixels with valid Combined calibration"); 259 MHCamera disp16(geomcam, "Saturation", "Pixels with saturated Hi Gain"); 260 MHCamera disp17(geomcam, "ConversionMeans", "Conversion HiGain.vs.LoGain Means"); 261 MHCamera disp18(geomcam, "ConversionSigmas", "Conversion HiGain.vs.LoGain Sigmas"); 262 MHCamera disp19(geomcam, "HiGainPickup", "Number Pickup events Hi Gain"); 263 MHCamera disp20(geomcam, "LoGainPickup", "Number Pickup events Lo Gain"); 264 MHCamera disp21(geomcam, "HiGainBlackout", "Number Blackout events Hi Gain"); 265 MHCamera disp22(geomcam, "LoGainBlackout", "Number Blackout events Lo Gain"); 266 MHCamera disp23(geomcam, "Excluded", "Pixels previously excluded"); 267 MHCamera disp24(geomcam, "UnSuitable", "Pixels not suited for further analysis"); 268 MHCamera disp25(geomcam, "UnReliable", "Pixels not reliable for further analysis"); 269 MHCamera disp26(geomcam, "HiGainOscillating", "Oscillating Pixels High Gain"); 270 MHCamera disp27(geomcam, "LoGainOscillating", "Oscillating Pixels Low Gain"); 271 MHCamera disp28(geomcam, "AbsTimeMean", "Abs. Arrival Times"); 272 MHCamera disp29(geomcam, "AbsTimeRms", "RMS of Arrival Times"); 273 MHCamera disp30(geomcam, "MeanTime", "Mean Rel. Arrival Times"); 274 MHCamera disp31(geomcam, "SigmaTime", "Sigma Rel. Arrival Times"); 275 MHCamera disp32(geomcam, "TimeProb", "Probability of Time Fit"); 276 MHCamera disp33(geomcam, "TimeNotFitValid", "Pixels with not valid fit results"); 277 MHCamera disp34(geomcam, "TimeOscillating", "Oscillating Pixels"); 304 278 305 279 // Fitted charge means and sigmas … … 665 639 Bool_t MJCalibration::FindColor() 666 640 { 667 668 const UInt_t nruns = fRuns->GetNumRuns(); 641 if (fSequence) 642 { 643 fColor = MCalibrationCam::kCT1; 644 return kTRUE; 645 } 646 647 const UInt_t nruns = fRuns->GetNumRuns(); 669 648 670 649 if (nruns == 0) … … 880 859 } 881 860 882 return kTRUE; 883 } 884 885 886 861 862 return kTRUE; 863 } 887 864 888 865 // -------------------------------------------------------------------------- … … 892 869 const char* MJCalibration::GetOutputFile() const 893 870 { 894 895 if (!fRuns)896 return "";897 898 return Form("%s/%s-F1.root", (const char*)fOutputPath, (const char*)fRuns->GetRunsAsFileName()); 899 } 900 871 if (fSequence) 872 return Form("%s/calib%06d.root", (const char*)fOutputPath, fSequence->GetSequence()); 873 if (!fRuns) 874 return ""; 875 876 return Form("%s/%s-F1.root", (const char*)fOutputPath, (const char*)fRuns->GetRunsAsFileName()); 877 } 901 878 902 879 Bool_t MJCalibration::IsUseBlindPixel() const 903 880 { 904 return TESTBIT(fDevices,kUseBlindPixel); 905 } 906 881 return TESTBIT(fDevices,kUseBlindPixel); 882 } 907 883 908 884 Bool_t MJCalibration::IsUsePINDiode() const 909 885 { 910 return TESTBIT(fDevices,kUsePINDiode); 911 } 912 913 914 915 886 return TESTBIT(fDevices,kUsePINDiode); 887 } 888 889 void MJCalibration::SetEnv(const char *env) 890 { 891 if (fEnv) 892 delete fEnv; 893 fEnv = new TEnv(env); 894 } 895 896 void MJCalibration::CheckEnv() 897 { 898 if (!fEnv) 899 return; 900 901 TString e1 = fEnv->GetValue("MJCalibration.OutputPath", ""); 902 if (!e1.IsNull()) 903 { 904 e1.ReplaceAll("\015", ""); 905 SetOutputPath(e1); 906 } 907 } 916 908 917 909 // -------------------------------------------------------------------------- … … 925 917 926 918 return kTRUE; 919 } 920 921 void MJCalibration::InitBlindPixel(MExtractBlindPixel &blindext, 922 MHCalibrationChargeBlindCam &blindcam) 923 { 924 Int_t run = fSequence ? fSequence->GetLastRun() : fRuns->GetRuns()[fRuns->GetNumRuns()-1]; 925 926 // 927 // Initialize the blind pixel. Unfortunately, there is a hardware difference 928 // in the first blind pixel until run "gkSecondBlindPixelInstallation" and the 929 // later setup. The first needs to use a filter because of the length of 930 // spurious NSB photon signals. The latter get better along extracting the amplitude 931 // from a small window. 932 // 933 if (run < gkSecondBlindPixelInstallation) 934 { 935 blindext.SetModified(kFALSE); 936 blindext.SetExtractionType(MExtractBlindPixel::kIntegral); 937 blindext.SetExtractionType(MExtractBlindPixel::kFilter); 938 blindext.SetRange(10,19,0,6); 939 blindext.SetNSBFilterLimit(70); 940 blindcam.SetFitFunc( MHCalibrationChargeBlindPix::kEPoisson5 ); 941 } 942 else 943 { 944 blindext.SetModified(kTRUE); 945 blindext.SetExtractionType(MExtractBlindPixel::kAmplitude); 946 blindext.SetExtractionType(MExtractBlindPixel::kFilter); 947 blindext.SetRange(5,8,0,2); 948 blindext.SetNSBFilterLimit(38); 949 950 if (run < gkThirdBlindPixelInstallation) 951 blindext.SetNumBlindPixels(2); 952 else 953 blindext.SetNumBlindPixels(3); 954 } 927 955 } 928 956 … … 965 993 Bool_t MJCalibration::ProcessFile(MPedestalCam &pedcam) 966 994 { 967 if (!fRuns) 968 { 969 *fLog << err << "No Runs choosen... abort." << endl; 970 return kFALSE; 971 } 995 if (!fRuns && !fSequence) 996 { 997 *fLog << err << "No Runs choosen... abort." << endl; 998 return kFALSE; 999 } 1000 1001 if (!fSequence && fRuns->GetNumRuns() != fRuns->GetNumEntries()) 1002 { 1003 *fLog << err << "Number of files found doesn't match number of runs... abort." 1004 << fRuns->GetNumRuns() << " vs. " << fRuns->GetNumEntries() << endl; 1005 return kFALSE; 1006 } 1007 1008 *fLog << inf; 1009 fLog->Separator(GetDescriptor()); 1010 1011 if (!FindColor()) 1012 return kFALSE; 1013 1014 *fLog << "Calculate MCalibrationCam from "; 1015 if (fSequence) 1016 *fLog << "Sequence #" << fSequence->GetSequence() << endl; 1017 else 1018 *fLog << "Runs " << fRuns->GetRunsAsString() << endl; 1019 *fLog << endl; 1020 1021 // Setup Tasklist 1022 MParList plist; 1023 MTaskList tlist; 1024 plist.AddToList(&tlist); 1025 1026 MReadMarsFile read("Events"); 1027 MRawFileRead rawread(NULL); 1028 1029 MDirIter iter; 1030 if (fSequence) 1031 fSequence->SetupCalRuns(iter); 1032 1033 if (fDataCheck) 1034 { 1035 rawread.AddFiles(fSequence ? iter : *fRuns); 1036 tlist.AddToList(&rawread); 1037 } 1038 else 1039 { 1040 read.DisableAutoScheme(); 1041 static_cast<MRead&>(read).AddFiles(fSequence ? iter : *fRuns); 1042 tlist.AddToList(&read); 1043 } 1044 1045 MHCalibrationChargeCam chargecam; 1046 MHCalibrationChargeBlindCam blindcam; 1047 1048 plist.AddToList(&pedcam); 1049 plist.AddToList(&chargecam); 1050 plist.AddToList(&blindcam); 1051 plist.AddToList(&fBadPixels); 1052 plist.AddToList(&fQECam); 1053 plist.AddToList(&fCalibrationCam); 1054 plist.AddToList(&fCalibrationBlindCam); 1055 plist.AddToList(&fCalibrationPINDiode); 1056 plist.AddToList(&fRelTimeCam); 1057 1058 MGeomApply apply; 1059 MBadPixelsMerge merge(&fBadPixels); 1060 MExtractPINDiode pinext; 1061 MExtractBlindPixel blindext; 1062 InitBlindPixel(blindext, blindcam); 1063 MExtractSlidingWindow extract2; 1064 MExtractTimeFastSpline timespline; 1065 MCalibrationChargeCalc calcalc; 1066 if (!fSequence) 1067 { 1068 calcalc.SetOutputPath(fOutputPath); 1069 calcalc.SetOutputFile(Form("%s-ChargeCalibStat.txt",(const char*)fRuns->GetRunsAsFileName())); 1070 } 1071 1072 if (fDebug) 1073 { 1074 chargecam.SetDebug(); 1075 calcalc.SetDebug(); 1076 } 1077 1078 // 1079 // As long as there are no DM's, have to colour by hand 1080 // 1081 calcalc.SetPulserColor(fColor); 1082 1083 MFillH fillpin("MHCalibrationChargePINDiode", "MExtractedSignalPINDiode"); 1084 MFillH fillbnd("MHCalibrationChargeBlindPix", "MExtractedSignalBlindPixel"); 1085 MFillH fillcam("MHCalibrationChargeCam", "MExtractedSignalCam"); 1086 MFillH filltme("MHCalibrationRelTimeCam", "MArrivalTimeCam"); 1087 fillpin.SetNameTab("PINDiode"); 1088 fillbnd.SetNameTab("BlindPix"); 1089 fillcam.SetNameTab("Charge"); 1090 filltme.SetNameTab("RelTimes"); 972 1091 973 if (fRuns->GetNumRuns() != fRuns->GetNumEntries()) 974 { 975 *fLog << err << "Number of files found doesn't match number of runs... abort." 976 << fRuns->GetNumRuns() << " vs. " << fRuns->GetNumEntries() << endl; 977 return kFALSE; 978 } 979 980 *fLog << inf; 981 fLog->Separator(GetDescriptor()); 982 983 if (!FindColor()) 984 return kFALSE; 985 986 *fLog << "Calculate MCalibrationCam from Runs " << fRuns->GetRunsAsString() << endl; 987 *fLog << endl; 988 989 // Setup Tasklist 990 MParList plist; 991 MTaskList tlist; 992 plist.AddToList(&tlist); 993 994 MReadMarsFile read("Events"); 995 MRawFileRead rawread(NULL); 996 997 if (fDataCheck) 998 { 999 rawread.AddFiles(*fRuns); 1000 tlist.AddToList(&rawread); 1001 } 1002 else 1003 { 1004 read.DisableAutoScheme(); 1005 static_cast<MRead&>(read).AddFiles(*fRuns); 1006 tlist.AddToList(&read); 1007 } 1008 1009 MHCalibrationChargeCam chargecam; 1010 MHCalibrationChargeBlindCam blindcam; 1011 1012 plist.AddToList(&pedcam); 1013 plist.AddToList(&chargecam); 1014 plist.AddToList(&blindcam); 1015 plist.AddToList(&fBadPixels); 1016 plist.AddToList(&fQECam); 1017 plist.AddToList(&fCalibrationCam); 1018 plist.AddToList(&fCalibrationBlindCam); 1019 plist.AddToList(&fCalibrationPINDiode); 1020 plist.AddToList(&fRelTimeCam); 1021 1022 MGeomApply apply; 1023 MExtractPINDiode pinext; 1024 MExtractBlindPixel blindext; 1025 1026 // 1027 // Initialize the blind pixel. Unfortunately, there is a hardware difference 1028 // in the first blind pixel until run "gkSecondBlindPixelInstallation" and the 1029 // later setup. The first needs to use a filter because of the length of 1030 // spurious NSB photon signals. The latter get better along extracting the amplitude 1031 // from a small window. 1032 // 1033 TArrayI arr = fRuns->GetRuns(); 1034 if (arr[fRuns->GetNumRuns()-1] < gkSecondBlindPixelInstallation) 1035 { 1036 blindext.SetModified(kFALSE); 1037 blindext.SetExtractionType(MExtractBlindPixel::kIntegral); 1038 blindext.SetExtractionType(MExtractBlindPixel::kFilter); 1039 blindext.SetRange(10,19,0,6); 1040 blindext.SetNSBFilterLimit(70); 1041 blindcam.SetFitFunc( MHCalibrationChargeBlindPix::kEPoisson5 ); 1042 } 1043 else 1044 { 1045 blindext.SetModified(kTRUE); 1046 blindext.SetExtractionType(MExtractBlindPixel::kAmplitude); 1047 blindext.SetExtractionType(MExtractBlindPixel::kFilter); 1048 blindext.SetRange(5,8,0,2); 1049 blindext.SetNSBFilterLimit(38); 1050 if (arr[fRuns->GetNumRuns()-1] < gkThirdBlindPixelInstallation) 1051 blindext.SetNumBlindPixels(2); 1052 else 1053 blindext.SetNumBlindPixels(3); 1054 } 1055 1056 MExtractSlidingWindow extract2; 1057 MExtractTimeFastSpline timespline; 1058 MCalibrationChargeCalc calcalc; 1059 calcalc.SetOutputPath(fOutputPath); 1060 calcalc.SetOutputFile(Form("%s-ChargeCalibStat.txt",(const char*)fRuns->GetRunsAsFileName())); 1061 1062 if (fDebug) 1063 { 1064 chargecam.SetDebug(); 1065 calcalc.SetDebug(); 1066 } 1067 1068 MCalibrationRelTimeCalc timecalc; 1069 timecalc.SetOutputPath(fOutputPath); 1070 timecalc.SetOutputFile(Form("%s-TimeCalibStat.txt",(const char*)fRuns->GetRunsAsFileName())); 1071 1072 // 1073 // As long as there are no DM's, have to colour by hand 1074 // 1075 calcalc.SetPulserColor(fColor); 1076 1077 MFillH fillpin("MHCalibrationChargePINDiode", "MExtractedSignalPINDiode"); 1078 MFillH fillbnd("MHCalibrationChargeBlindCam", "MExtractedSignalBlindPixel"); 1079 MFillH fillcam("MHCalibrationChargeCam", "MExtractedSignalCam"); 1080 MFillH filltme("MHCalibrationRelTimeCam", "MArrivalTimeCam"); 1081 fillpin.SetNameTab("PINDiode"); 1082 fillbnd.SetNameTab("BlindPix"); 1083 fillcam.SetNameTab("Charge"); 1084 filltme.SetNameTab("RelTimes"); 1085 1086 TString drawoption; 1087 1088 if (fDisplayType == kDataCheckDisplay) 1089 drawoption += "datacheck"; 1090 if (fDisplayType == kFullDisplay) 1091 drawoption += " all"; 1092 1093 fillcam.SetDrawOption(drawoption.Data()); 1094 fillbnd.SetDrawOption(drawoption.Data()); 1095 fillpin.SetDrawOption(drawoption.Data()); 1096 filltme.SetDrawOption(drawoption.Data()); 1097 // 1098 // Apply a filter against cosmics 1099 // (will have to be needed in the future 1100 // when the calibration hardware-trigger is working) 1101 // 1102 MFCosmics cosmics; 1103 MContinue cont(&cosmics); 1104 1105 // tlist.AddToList(&merge); 1106 tlist.AddToList(&apply); 1107 1108 if (fExtractor) 1109 tlist.AddToList(fExtractor); 1110 else 1111 { 1112 *fLog << warn << GetDescriptor() 1113 << ": No extractor has been chosen, take default MExtractSlidingWindow " << endl; 1114 tlist.AddToList(&extract2); 1115 } 1116 1117 tlist.AddToList(&pinext); 1118 tlist.AddToList(&blindext); 1119 1120 if (fRelTimes) 1121 { 1122 if (fTimeExtractor) 1123 tlist.AddToList(fTimeExtractor); 1124 else 1092 TString drawoption; 1093 1094 if (fDisplayType == kDataCheckDisplay) 1095 drawoption += "datacheck"; 1096 if (fDisplayType == kFullDisplay) 1097 drawoption += " all"; 1098 1099 fillcam.SetDrawOption(drawoption.Data()); 1100 fillbnd.SetDrawOption(drawoption.Data()); 1101 fillpin.SetDrawOption(drawoption.Data()); 1102 filltme.SetDrawOption(drawoption.Data()); 1103 1104 // 1105 // Apply a filter against cosmics 1106 // (will have to be needed in the future 1107 // when the calibration hardware-trigger is working) 1108 // 1109 MFCosmics cosmics; 1110 MContinue cont(&cosmics); 1111 1112 tlist.AddToList(&merge); 1113 tlist.AddToList(&apply); 1114 1115 MTaskEnv taskenv("ExtractSignal"); 1116 taskenv.SetDefault(fExtractor ? fExtractor : &extract2); 1117 1118 tlist.AddToList(&taskenv); 1119 tlist.AddToList(&pinext); 1120 tlist.AddToList(&blindext); 1121 1122 MTaskEnv taskenv2("ExtractTime"); 1123 taskenv2.SetDefault(fTimeExtractor ? fTimeExtractor : ×pline); 1124 1125 if (fRelTimes) 1126 tlist.AddToList(&taskenv2); 1127 1128 if (fColor == MCalibrationCam::kCT1) 1129 tlist.AddToList(&cont); 1130 1131 tlist.AddToList(&fillcam); 1132 tlist.AddToList(&calcalc); 1133 1134 MCalibrationRelTimeCalc timecalc; 1135 if (fRelTimes) 1136 { 1137 tlist.AddToList(&filltme); 1138 tlist.AddToList(&timecalc); 1139 } 1140 1141 1142 // Create and setup the eventloop 1143 MEvtLoop evtloop(fName); 1144 evtloop.SetParList(&plist); 1145 evtloop.SetDisplay(fDisplay); 1146 evtloop.SetLogStream(fLog); 1147 if (fEnv) 1148 evtloop.ReadEnv(*fEnv); 1149 1150 // Execute first analysis 1151 if (!evtloop.Eventloop()) 1152 { 1153 *fLog << err << GetDescriptor() << ": Failed." << endl; 1154 return kFALSE; 1155 } 1156 1157 tlist.PrintStatistics(); 1158 1159 // 1160 // The next lines are necessary in order to avoid that 1161 // the last entry drawn by MFillH gets deleted again from 1162 // the display. No idea where this comes from... 1163 // 1164 /* 1165 if (fDisplay) 1166 { 1167 if (IsUsePINDiode()) 1125 1168 { 1126 *fLog << warn << GetDescriptor()1127 << ": No extractor has been chosen, take default MTimeExtractSpline " << endl;1128 tlist.AddToList(×pline);1169 MHCalibrationChargePINDiode *pin = 1170 (MHCalibrationChargePINDiode*)plist.FindObject("MHCalibrationChargePINDiode"); 1171 pin->DrawClone(Form("nonew %s",drawoption.Data())); 1129 1172 } 1130 } 1131 1132 if (fColor == MCalibrationCam::kCT1) 1133 tlist.AddToList(&cont); 1134 1135 tlist.AddToList(&fillcam); 1136 1137 if (fRelTimes) 1138 { 1139 tlist.AddToList(&filltme); 1140 tlist.AddToList(&timecalc); 1141 } 1142 1143 if (IsUseBlindPixel()) 1144 tlist.AddToList(&fillbnd); 1145 if (IsUsePINDiode()) 1146 tlist.AddToList(&fillpin); 1147 1148 tlist.AddToList(&calcalc); 1149 1150 // Create and setup the eventloop 1151 MEvtLoop evtloop(fName); 1152 evtloop.SetParList(&plist); 1153 evtloop.SetDisplay(fDisplay); 1154 evtloop.SetLogStream(fLog); 1155 1156 // Execute first analysis 1157 if (!evtloop.Eventloop()) 1158 { 1159 *fLog << err << GetDescriptor() << ": Failed." << endl; 1160 return kFALSE; 1161 } 1162 1163 tlist.PrintStatistics(); 1164 1165 // 1166 // The next lines are necessary in order to avoid that 1167 // the last entry drawn by MFillH gets deleted again from 1168 // the display. No idea where this comes from... 1169 // 1170 if (fDisplay) 1171 { 1172 if (IsUsePINDiode()) 1173 else if (IsUseBlindPixel()) 1173 1174 { 1174 MHCalibrationChargePINDiode *pin =1175 (MHCalibrationChargePINDiode*)plist.FindObject("MHCalibrationChargePINDiode");1176 pin->DrawClone(Form("nonew %s",drawoption.Data()));1175 MHCalibrationChargeBlindCam *cam = 1176 (MHCalibrationChargeBlindCam*)plist.FindObject("MHCalibrationChargeBlindCam"); 1177 cam->DrawClone(Form("nonew %s",drawoption.Data())); 1177 1178 } 1178 else if (IsUseBlindPixel())1179 else if (fRelTimes) 1179 1180 { 1180 MHCalibrationChargeBlindCam *cam =1181 (MHCalibrationChargeBlindCam*)plist.FindObject("MHCalibrationChargeBlindCam");1182 cam->DrawClone(Form("nonew %s",drawoption.Data()));1181 MHCalibrationRelTimeCam *cam = 1182 (MHCalibrationRelTimeCam*)plist.FindObject("MHCalibrationRelTimeCam"); 1183 cam->DrawClone(Form("nonew %s",drawoption.Data())); 1183 1184 } 1184 else if (fRelTimes)1185 else 1185 1186 { 1186 MHCalibrationRelTimeCam *cam =1187 (MHCalibrationRelTimeCam*)plist.FindObject("MHCalibrationRelTimeCam");1188 cam->DrawClone(Form("nonew %s",drawoption.Data()));1187 MHCalibrationChargeCam *cam = 1188 (MHCalibrationChargeCam*)plist.FindObject("MHCalibrationChargeCam"); 1189 cam->DrawClone(Form("nonew %s",drawoption.Data())); 1189 1190 } 1190 else 1191 { 1192 MHCalibrationChargeCam *cam = 1193 (MHCalibrationChargeCam*)plist.FindObject("MHCalibrationChargeCam"); 1194 cam->DrawClone(Form("nonew %s",drawoption.Data())); 1195 } 1196 } 1197 1198 DisplayResult(plist); 1199 1200 if (!WriteResult()) 1201 return kFALSE; 1202 1203 *fLog << inf << GetDescriptor() << ": Done." << endl; 1204 1205 return kTRUE; 1191 } 1192 */ 1193 1194 DisplayResult(plist); 1195 1196 if (!WriteResult()) 1197 return kFALSE; 1198 1199 *fLog << inf << GetDescriptor() << ": Done." << endl; 1200 1201 return kTRUE; 1206 1202 } 1207 1203 -
trunk/MagicSoft/Mars/mjobs/MJCalibration.h
r4604 r4609 21 21 #endif 22 22 23 class TEnv; 23 24 class MRunIter; 25 class MSequence; 24 26 class MParList; 25 27 class MPedestalCam; 26 28 class MExtractor; 27 29 class MExtractTime; 30 31 class MExtractBlindPixel; 32 class MHCalibrationChargeBlindCam; 33 28 34 class MJCalibration : public MParContainer 29 35 { 30 36 private: 37 static const Int_t gkIFAEBoxInaugurationRun; // Run number of first IFAE box calibration 38 static const Int_t gkSecondBlindPixelInstallation; // Run number upon which second blind pixel was installed 39 static const Int_t gkThirdBlindPixelInstallation; // Run number upon which third blind pixel was installed 31 40 32 static const Int_t gkIFAEBoxInaugurationRun; // Run number of first IFAE box calibration 33 static const Int_t gkSecondBlindPixelInstallation; // Run number upon which second blind pixel was installed 34 static const Int_t gkThirdBlindPixelInstallation; // Run number upon which third blind pixel was installed 41 TString fOutputPath; // Path to the output files 35 42 36 TString fOutputPath; // Path to the output files 43 TEnv *fEnv; 44 MRunIter *fRuns; // Calibration files 45 MSequence *fSequence; // Sequence 46 47 MExtractor *fExtractor; // Signal extractor 48 MExtractTime *fTimeExtractor; // Arrival Time extractor 49 50 MBadPixelsCam fBadPixels; // Bad Pixels cam, can be set from previous runs 51 MCalibrationChargeCam fCalibrationCam; // Calibration conversion factors FADC2Phe 52 MCalibrationChargeBlindCam fCalibrationBlindCam; // Calibration from Blind Pixel(s) 53 MCalibrationChargePINDiode fCalibrationPINDiode; // Calibration from PIN Diode 54 MCalibrationQECam fQECam; // Quantum efficiency, can be set from previous runs 55 MCalibrationRelTimeCam fRelTimeCam; // Calibration constants rel. times 56 57 MCalibrationCam::PulserColor_t fColor; // Colour of the pulsed LEDs 58 59 enum Display_t // Possible Display types 60 { 61 kFullDisplay, 62 kDataCheckDisplay, 63 kNormalDisplay 64 }; 65 66 Display_t fDisplayType; // Chosen Display type 67 68 enum Device_t // Possible devices for calibration 69 { 70 kUseBlindPixel, 71 kUsePINDiode 72 }; 73 74 Byte_t fDevices; // Bit-field for used devices for calibration 75 76 Bool_t fRelTimes; // Flag if relative times have to be calibrated 77 Bool_t fDataCheck; // Flag if the data check is run on raw data 78 Bool_t fDebug; 79 80 void DisplayResult(MParList &plist); 81 Bool_t WriteResult(); 82 void CheckEnv(); 83 84 // WORKAROUNDS!!! 85 Bool_t FindColor(); 86 void InitBlindPixel(MExtractBlindPixel &blindext, 87 MHCalibrationChargeBlindCam &blindcam); 88 89 public: 90 MJCalibration(const char *name=NULL, const char *title=NULL); 91 ~MJCalibration(); 37 92 38 MRunIter *fRuns; // Calibration files 39 MExtractor *fExtractor; // Signal extractor 40 MExtractTime *fTimeExtractor; // Arrival Time extractor 93 const char* GetOutputFile() const; 94 void SetEnv(const char *env); 95 96 MCalibrationChargeCam &GetCalibrationCam() { return fCalibrationCam; } 97 MCalibrationRelTimeCam &GetRelTimeCam() { return fRelTimeCam; } 98 MCalibrationQECam &GetQECam() { return fQECam; } 99 MBadPixelsCam &GetBadPixels() { return fBadPixels; } 100 101 Bool_t IsUseBlindPixel() const; 102 Bool_t IsUsePINDiode() const; 103 104 void SetBadPixels(const MBadPixelsCam &bad) { bad.Copy(fBadPixels); } 105 void SetExtractor(MExtractor* ext) { fExtractor = ext; } 106 void SetTimeExtractor(MExtractTime* ext) { fTimeExtractor = ext; } 107 void SetQECam(const MCalibrationQECam &qe) { qe.Copy(fQECam); } 108 void SetColor(const MCalibrationCam::PulserColor_t color) { fColor = color; } 109 110 void SetInput(MRunIter *iter) { fRuns=iter; } 111 void SetSequence(MSequence *seq) { fSequence=seq; } 112 void SetOutputPath(const char *path="."); 113 114 // Displays 115 void SetFullDisplay() { fDisplayType = kFullDisplay; } 116 void SetDataCheckDisplay() { fDisplayType = kDataCheckDisplay; } 117 void SetNormalDisplay() { fDisplayType = kNormalDisplay; } 118 119 // Rel. Time 120 void SetRelTimeCalibration(const Bool_t b=kTRUE) { fRelTimes = b; } 121 122 // Data Check 123 void SetDataCheck(const Bool_t b=kTRUE) { fDataCheck = b; SetDataCheckDisplay(); } 124 125 // Debug 126 void SetDebug(const Bool_t b=kTRUE) { fDebug = b; } 127 128 // Devices 129 void SetUseBlindPixel(const Bool_t b=kTRUE); 130 void SetUsePINDiode(const Bool_t b=kTRUE); 41 131 42 MBadPixelsCam fBadPixels; // Bad Pixels cam, can be set from previous runs 43 MCalibrationChargeCam fCalibrationCam; // Calibration conversion factors FADC2Phe 44 MCalibrationChargeBlindCam fCalibrationBlindCam; // Calibration from Blind Pixel(s) 45 MCalibrationChargePINDiode fCalibrationPINDiode; // Calibration from PIN Diode 46 MCalibrationQECam fQECam; // Quantum efficiency, can be set from previous runs 47 MCalibrationRelTimeCam fRelTimeCam; // Calibration constants rel. times 132 Bool_t ReadCalibrationCam(); 133 Bool_t ProcessFile(MPedestalCam &pedcam); 134 Bool_t Process(MPedestalCam &pedcam); 48 135 49 MCalibrationCam::PulserColor_t fColor; // Colour of the pulsed LEDs 50 51 enum Display_t { kFullDisplay, kDataCheckDisplay, kNormalDisplay }; // Possible Display types 52 53 Display_t fDisplayType; // Chosen Display type 54 55 enum Device_t { kUseBlindPixel, kUsePINDiode }; // Possible devices for calibration 56 57 Byte_t fDevices; // Bit-field for used devices for calibration 58 59 Bool_t fRelTimes; // Flag if relative times have to be calibrated 60 Bool_t fDataCheck; // Flag if the data check is run on raw data 61 Bool_t fDebug; 62 63 void DisplayResult(MParList &plist); 64 Bool_t WriteResult(); 65 Bool_t FindColor(); 66 67 public: 68 69 MJCalibration(const char *name=NULL, const char *title=NULL); 70 ~MJCalibration() {} 71 72 const char* GetOutputFile() const; 73 74 MCalibrationChargeCam &GetCalibrationCam() { return fCalibrationCam; } 75 MCalibrationRelTimeCam &GetRelTimeCam() { return fRelTimeCam; } 76 MCalibrationQECam &GetQECam() { return fQECam; } 77 MBadPixelsCam &GetBadPixels() { return fBadPixels; } 78 79 Bool_t IsUseBlindPixel() const; 80 Bool_t IsUsePINDiode() const; 81 82 void SetBadPixels(const MBadPixelsCam &bad) { bad.Copy(fBadPixels); } 83 void SetExtractor(MExtractor* ext) { fExtractor = ext; } 84 void SetTimeExtractor(MExtractTime* ext) { fTimeExtractor = ext; } 85 void SetQECam (const MCalibrationQECam &qe) { qe.Copy(fQECam); } 86 void SetColor (const MCalibrationCam::PulserColor_t color) { fColor = color; } 87 88 void SetInput(MRunIter *iter) { fRuns=iter; } 89 void SetOutputPath(const char *path="."); 90 91 // Displays 92 void SetFullDisplay() { fDisplayType = kFullDisplay; } 93 void SetDataCheckDisplay() { fDisplayType = kDataCheckDisplay; } 94 void SetNormalDisplay() { fDisplayType = kNormalDisplay; } 95 96 // Rel. Time 97 void SetRelTimeCalibration(const Bool_t b=kTRUE) { fRelTimes = b; } 98 99 // Data Check 100 void SetDataCheck (const Bool_t b=kTRUE) { fDataCheck = b; SetDataCheckDisplay(); } 101 102 // Debug 103 void SetDebug (const Bool_t b=kTRUE) { fDebug = b; } 104 105 // Devices 106 void SetUseBlindPixel( const Bool_t b=kTRUE ); 107 void SetUsePINDiode ( const Bool_t b=kTRUE ); 108 109 Bool_t ReadCalibrationCam(); 110 Bool_t ProcessFile( MPedestalCam &pedcam ); 111 Bool_t Process ( MPedestalCam &pedcam ); 112 113 ClassDef(MJCalibration, 0) // Tool to run a calibration per pulser colour and intensity 136 ClassDef(MJCalibration, 0) // Tool to run a calibration per pulser colour and intensity 114 137 }; 115 138 -
trunk/MagicSoft/Mars/mjobs/MJPedestal.cc
r4601 r4609 104 104 { 105 105 if (fSequence) 106 return Form("%s/calped%0 5d", (const char*)fOutputPath, fSequence->GetSequence());106 return Form("%s/calped%06d.root", (const char*)fOutputPath, fSequence->GetSequence()); 107 107 108 108 if (!fRuns) … … 422 422 423 423 return kTRUE; 424 } 425 426 void MJPedestal::CheckEnv() 427 { 428 if (!fEnv) 429 return; 430 431 TString e1 = fEnv->GetValue("MJPedestal.OutputPath", ""); 432 if (!e1.IsNull()) 433 { 434 e1.ReplaceAll("\015", ""); 435 SetOutputPath(e1); 436 } 424 437 } 425 438 … … 449 462 *fLog << endl; 450 463 464 CheckEnv(); 465 451 466 MParList plist; 452 467 MTaskList tlist; … … 493 508 pedcalc.SetRange(fExtractor->GetHiGainFirst(), fExtractor->GetHiGainLast()); 494 509 } 510 /* 495 511 else 496 512 { … … 498 514 *fLog << ": No extractor has been chosen, take default number of FADC samples " << endl; 499 515 } 516 */ 500 517 501 518 tlist.AddToList(&geomapl); -
trunk/MagicSoft/Mars/mjobs/MJPedestal.h
r4601 r4609 58 58 void DisplayOutliers(TH1D *hist) const; 59 59 60 void CheckEnv(); 61 60 62 public: 61 63 MJPedestal(const char *name=NULL, const char *title=NULL); -
trunk/MagicSoft/Mars/mjobs/MSequence.h
r4601 r4609 50 50 // Getter 51 51 UInt_t GetSequence() const { return fSequence; } 52 UInt_t GetLastRun() const { return fLastRun; } 52 53 53 54 ClassDef(MSequence, 0) -
trunk/MagicSoft/Mars/mpedestal/MPedCalcFromLoGain.cc
r4601 r4609 185 185 fGeom(NULL), fPedContainerName("MPedestalCam") 186 186 { 187 fName = name ? name : "MPedCalcFromLoGain"; 188 fTitle = title ? title : "Task to calculate pedestals from pedestal runs raw data"; 189 190 AddToBranchList("fHiGainPixId"); 191 AddToBranchList("fHiGainFadcSamples"); 192 193 SetRange(fgHiGainFirst, fgHiGainLast, fgLoGainFirst, fgLoGainLast); 194 195 SetMaxHiGainVar(fgMaxHiGainVar); 196 SetPedestalUpdate(kTRUE); 197 Clear(); 187 fName = name ? name : "MPedCalcFromLoGain"; 188 fTitle = title ? title : "Task to calculate pedestals from pedestal runs raw data"; 189 190 AddToBranchList("fHiGainPixId"); 191 AddToBranchList("fLoGainPixId"); 192 AddToBranchList("fHiGainFadcSamples"); 193 AddToBranchList("fLoGainFadcSamples"); 194 195 SetRange(); 196 SetMaxHiGainVar(); 197 SetPedestalUpdate(kTRUE); 198 199 Clear(); 198 200 } 199 201 … … 208 210 void MPedCalcFromLoGain::Clear(const Option_t *o) 209 211 { 210 fRawEvt = NULL;211 fRunHeader = NULL;212 fPedestals = NULL;213 214 // If the size is yet set, set the size215 if (fSumx.GetSize()>0)216 {217 // Reset contents of arrays.218 fSumx.Reset();219 fSumx2.Reset();220 fSumAB0.Reset();221 fSumAB1.Reset();222 fNumEventsUsed.Reset();223 fTotalCounter.Reset();224 }212 fRawEvt = NULL; 213 fRunHeader = NULL; 214 fPedestals = NULL; 215 216 // If the size is yet set, set the size 217 if (fSumx.GetSize()>0) 218 { 219 // Reset contents of arrays. 220 fSumx.Reset(); 221 fSumx2.Reset(); 222 fSumAB0.Reset(); 223 fSumAB1.Reset(); 224 fNumEventsUsed.Reset(); 225 fTotalCounter.Reset(); 226 } 225 227 } 226 228 … … 235 237 void MPedCalcFromLoGain::SetRange(Byte_t hifirst, Byte_t hilast, Byte_t lofirst, Byte_t lolast) 236 238 { 237 MExtractor::SetRange(hifirst,hilast,lofirst,lolast);238 239 //240 // Redo the checks if the window is still inside the ranges241 //242 SetWindowSize(fWindowSizeHiGain,fWindowSizeLoGain);239 MExtractor::SetRange(hifirst, hilast, lofirst, lolast); 240 241 // 242 // Redo the checks if the window is still inside the ranges 243 // 244 SetWindowSize(fWindowSizeHiGain,fWindowSizeLoGain); 243 245 } 244 246 … … 247 249 void MPedCalcFromLoGain::SetMaxHiGainVar(Byte_t maxvar) 248 250 { 249 fMaxHiGainVar = maxvar;251 fMaxHiGainVar = maxvar; 250 252 } 251 253 … … 265 267 void MPedCalcFromLoGain::SetWindowSize(Byte_t windowh, Byte_t windowl) 266 268 { 267 fWindowSizeHiGain = windowh & ~1; 268 fWindowSizeLoGain = windowl & ~1; 269 270 if (fWindowSizeHiGain != windowh) 271 *fLog << warn << GetDescriptor() << ": Hi Gain window size has to be even, set to: " 272 << int(fWindowSizeHiGain) << " samples " << endl; 269 fWindowSizeHiGain = windowh & ~1; 270 fWindowSizeLoGain = windowl & ~1; 273 271 274 if (fWindowSizeLoGain != windowl) 275 *fLog << warn << GetDescriptor() << ": Lo Gain window size has to be even, set to: " 276 << int(fWindowSizeLoGain) << " samples " << endl; 272 if (fWindowSizeHiGain != windowh) 273 { 274 *fLog << warn; 275 *fLog << GetDescriptor() << ": HiGain window has to be even, set to: "; 276 *fLog << int(fWindowSizeHiGain) << " samples " << endl; 277 } 277 278 278 const Byte_t availhirange = (fHiGainLast-fHiGainFirst+1) & ~1; 279 const Byte_t availlorange = (fLoGainLast-fLoGainFirst+1) & ~1; 280 281 if (fWindowSizeHiGain > availhirange) 282 { 283 *fLog << warn << GetDescriptor() 284 << Form(": Hi Gain window size: %2i is bigger than available range: [%2i,%2i]", 285 (int)fWindowSizeHiGain, (int)fHiGainFirst, (int)fHiGainLast) << endl; 286 *fLog << warn << GetDescriptor() 287 << ": Will set window size to: " << (int)availhirange << endl; 288 fWindowSizeHiGain = availhirange; 289 } 279 if (fWindowSizeLoGain != windowl) 280 { 281 *fLog << warn; 282 *fLog << GetDescriptor() << ": Lo Gain window has to be even, set to: "; 283 *fLog << int(fWindowSizeLoGain) << " samples " << endl; 284 } 285 286 const Byte_t availhirange = (fHiGainLast-fHiGainFirst+1) & ~1; 287 const Byte_t availlorange = (fLoGainLast-fLoGainFirst+1) & ~1; 290 288 291 if (fWindowSizeLoGain > availlorange) 292 { 293 *fLog << warn << GetDescriptor() 294 << Form(": Lo Gain window size: %2i is bigger than available range: [%2i,%2i]", 295 (int)fWindowSizeLoGain, (int)fLoGainFirst, (int)fLoGainLast) << endl; 296 *fLog << warn << GetDescriptor() 297 << ": Will set window size to: " << (int)availlorange << endl; 298 fWindowSizeLoGain = availlorange; 299 } 300 301 fNumHiGainSamples = (Float_t)fWindowSizeHiGain; 302 fNumLoGainSamples = (Float_t)fWindowSizeLoGain; 303 304 fSqrtHiGainSamples = TMath::Sqrt(fNumHiGainSamples); 305 fSqrtLoGainSamples = TMath::Sqrt(fNumLoGainSamples); 289 if (fWindowSizeHiGain > availhirange) 290 { 291 *fLog << warn; 292 *fLog << GetDescriptor() << ": HiGain window " << (int)fWindowSizeHiGain; 293 *fLog << " out of range [" << (int)fHiGainFirst; 294 *fLog << "," << (int)fHiGainLast << "]" << endl; 295 *fLog << "Will set window size to " << (int)availhirange << endl; 296 fWindowSizeHiGain = availhirange; 297 } 298 299 if (fWindowSizeLoGain > availlorange) 300 { 301 *fLog << warn; 302 *fLog << GetDescriptor() << ": LoGain window " << (int)fWindowSizeLoGain; 303 *fLog << " out of range [" << (int)fLoGainFirst; 304 *fLog << "," << (int)fLoGainLast << "]" << endl; 305 *fLog << "Will set window size to " << (int)availlorange << endl; 306 fWindowSizeLoGain = availlorange; 307 } 308 /* 309 fNumHiGainSamples = (Float_t)fWindowSizeHiGain; 310 fNumLoGainSamples = (Float_t)fWindowSizeLoGain; 311 312 fSqrtHiGainSamples = TMath::Sqrt(fNumHiGainSamples); 313 fSqrtLoGainSamples = TMath::Sqrt(fNumLoGainSamples); 314 */ 306 315 } 307 316 … … 369 378 Bool_t MPedCalcFromLoGain::ReInit(MParList *pList) 370 379 { 371 Int_t lastdesired = (Int_t)fLoGainLast; 372 Int_t lastavailable = (Int_t)fRunHeader->GetNumSamplesLoGain()-1; 373 374 if (lastdesired > lastavailable) 375 { 376 const Int_t diff = lastdesired - lastavailable; 377 *fLog << endl; 378 *fLog << warn << GetDescriptor() 379 << Form(": Selected Lo Gain FADC Window [%2i,%2i] ranges out of the available limits: [0,%2i].", 380 (int)fLoGainFirst, lastdesired, lastavailable) << endl; 381 *fLog << GetDescriptor() << ": Will reduce the upper edge to " << (int)(fLoGainLast - diff) << endl; 382 SetRange(fHiGainFirst, fHiGainLast, fLoGainFirst, fLoGainLast-diff); 383 } 384 385 lastdesired = (Int_t)fHiGainLast; 386 lastavailable = (Int_t)fRunHeader->GetNumSamplesHiGain()-1; 387 388 if (lastdesired > lastavailable) 389 { 390 const Int_t diff = lastdesired - lastavailable; 391 *fLog << endl; 392 *fLog << warn << GetDescriptor() 393 << Form(": Selected Hi Gain Range [%2i,%2i] ranges out of the available limits: [0,%2i].", 394 (int)fHiGainFirst, lastdesired, lastavailable) << endl; 395 *fLog << warn << GetDescriptor() 396 << Form(": Will possibly use %2i samples from the Low-Gain for the High-Gain range", diff) 397 << endl; 398 fHiGainLast -= diff; 399 fHiLoLast = diff; 400 } 401 402 lastdesired = (Int_t)fHiGainFirst+fWindowSizeHiGain-1; 403 lastavailable = (Int_t)fRunHeader->GetNumSamplesHiGain()-1; 404 405 if (lastdesired > lastavailable) 406 { 407 const Int_t diff = lastdesired - lastavailable; 408 *fLog << endl; 409 *fLog << warn << GetDescriptor() 410 << Form(": Selected Hi Gain FADC Window size %2i ranges out of the available limits: [0,%2i].", 411 (int)fWindowSizeHiGain, lastavailable) << endl; 412 *fLog << warn << GetDescriptor() 413 << Form(": Will use %2i samples from the Low-Gain for the High-Gain extraction", diff) 414 << endl; 415 416 if ((Int_t)fWindowSizeHiGain > diff) 417 { 418 fWindowSizeHiGain -= diff; 419 fWindowSizeLoGain += diff; 420 } 421 else 422 { 423 fWindowSizeLoGain += fWindowSizeHiGain; 424 fLoGainFirst = diff-fWindowSizeHiGain; 425 fWindowSizeHiGain = 0; 426 } 427 } 428 429 430 // If the size is not yet set, set the size 431 if (fSumx.GetSize()==0) 432 { 433 const Int_t npixels = fPedestals->GetSize(); 434 435 fSumx. Set(npixels); 436 fSumx2.Set(npixels); 437 fSumAB0.Set(npixels); 438 fSumAB1.Set(npixels); 439 fNumEventsUsed.Set(npixels); 440 fTotalCounter.Set(npixels); 441 442 // Reset contents of arrays. 443 fSumx.Reset(); 444 fSumx2.Reset(); 445 fSumAB0.Reset(); 446 fSumAB1.Reset(); 447 fNumEventsUsed.Reset(); 448 fTotalCounter.Reset(); 449 } 450 451 if (fWindowSizeHiGain==0 && fWindowSizeLoGain==0) 452 { 453 *fLog << err << GetDescriptor() 454 << ": Number of extracted Slices is 0, abort ... " << endl; 455 return kFALSE; 456 } 457 458 *fLog << endl; 459 *fLog << inf << GetDescriptor() << ": Taking " << (int)fWindowSizeHiGain 460 << " HiGain FADC samples starting with slice: " << (int)fHiGainFirst << endl; 461 *fLog << inf << GetDescriptor() << ": Taking " << (int)fWindowSizeLoGain 462 << " LoGain FADC samples starting with slice: " << (int)fLoGainFirst << endl; 463 464 return kTRUE; 380 if (fRunHeader->GetNumSamplesHiGain()<1) 381 { 382 *fLog << err << "ERROR - Number of available HiGainSamples<1... abort." << endl; 383 return kFALSE; 384 } 385 if (fRunHeader->GetNumSamplesLoGain()<1) 386 { 387 *fLog << err << "ERROR - Number of available LoGainSamples<1... abort." << endl; 388 return kFALSE; 389 } 390 391 fHiGainFirst = TMath::Max(fHiGainFirst, 0); 392 fHiGainLast = TMath::Min(fHiGainLast, fRunHeader->GetNumSamplesHiGain()-1); 393 394 fLoGainFirst = TMath::Max(fLoGainFirst, 0); 395 fLoGainLast = TMath::Min(fLoGainLast, fRunHeader->GetNumSamplesLoGain()-1); 396 397 fWindowSizeHiGain = TMath::Min(fWindowSizeHiGain, fHiGainLast-fHiGainFirst+1); 398 fWindowSizeLoGain = TMath::Min(fWindowSizeLoGain, fLoGainLast-fLoGainFirst+1); 399 400 SetRange(fHiGainFirst, fHiGainLast, fLoGainFirst, fLoGainLast); 401 402 *fLog << inf << endl; 403 *fLog << "Taking " << Form("%2d", (int)fWindowSizeHiGain) << " HiGain from " << (int)fHiGainFirst << endl; 404 *fLog << "Taking " << Form("%2d", (int)fWindowSizeLoGain) << " LoGain from " << (int)fLoGainFirst << endl; 405 406 if (fWindowSizeHiGain==0) 407 { 408 *fLog << err << "ERROR - HiGain windows size == 0... abort." << endl; 409 return kFALSE; 410 } 411 if (fWindowSizeLoGain==0) 412 { 413 *fLog << err << "ERROR - HiGain windows size == 0... abort." << endl; 414 return kFALSE; 415 } 416 417 // If the size is not yet set, set the size 418 if (fSumx.GetSize()==0) 419 { 420 const Int_t npixels = fPedestals->GetSize(); 421 422 fSumx. Set(npixels); 423 fSumx2.Set(npixels); 424 fSumAB0.Set(npixels); 425 fSumAB1.Set(npixels); 426 fNumEventsUsed.Set(npixels); 427 fTotalCounter.Set(npixels); 428 429 // Reset contents of arrays. 430 fSumx.Reset(); 431 fSumx2.Reset(); 432 fSumAB0.Reset(); 433 fSumAB1.Reset(); 434 fNumEventsUsed.Reset(); 435 fTotalCounter.Reset(); 436 } 437 438 return kTRUE; 465 439 } 466 440 … … 515 489 // Find the maximum and minimum signal per slice in the high gain window 516 490 do { 517 if (*ptr > max) {491 if (*ptr > max) 518 492 max = *ptr; 519 } 520 if (*ptr < min) { 493 if (*ptr < min) 521 494 min = *ptr; 522 }523 495 } while (++ptr != end); 524 496 … … 527 499 continue; 528 500 529 Byte_t *firstSlice = pixel.GetLoGainSamples() + fLoGainFirst; 530 Byte_t *lastSlice = firstSlice + fWindowSizeLoGain; 531 532 Byte_t *slice = firstSlice; 501 ptr = pixel.GetLoGainSamples() + fLoGainFirst; 502 end = ptr + fWindowSizeLoGain; 503 504 Byte_t *firstSlice = ptr; 505 533 506 do { 534 sum += * slice;535 sqr += * slice * *slice;536 } while (++ slice != lastSlice);507 sum += *ptr; 508 sqr += *ptr * *ptr; 509 } while (++ptr != end); 537 510 538 511 const Float_t msum = (Float_t)sum; … … 569 542 570 543 if (fPedestalUpdate) 544 { 545 fPedestals->ReCalc(*fGeom); 571 546 fPedestals->SetReadyToSave(); 547 } 572 548 573 549 return kTRUE; … … 581 557 { 582 558 // Compute pedestals and rms from the whole run 583 if (fPedestalUpdate )559 if (fPedestalUpdate || GetNumExecutions()<1) 584 560 return kTRUE; 585 561 586 *fLog << flush << inf << "Calculating Pedestals..." << flush; 562 *fLog << flush << inf << "Calculating pedestals..." << flush; 563 564 Double_t sum = 0; 587 565 588 566 const Int_t npix = fGeom->GetNumPixels(); … … 592 570 if (n>1) 593 571 Calc(n, idx); 594 } 595 572 sum += n; 573 } 574 575 *fLog << flush << inf << "Calculating means..." << flush; 576 577 fPedestals->SetTotalEntries((UInt_t)(sum/npix*(fWindowSizeLoGain+fWindowSizeHiGain))); 578 fPedestals->ReCalc(*fGeom); 596 579 fPedestals->SetReadyToSave(); 597 580 return kTRUE; … … 621 604 SetWindowSize(hw, lw); 622 605 623 Int_t num = fNumEventsDump;624 606 if (IsEnvDefined(env, prefix, "NumEventsDump", print)) 625 607 { 626 num = GetEnvValue(env, prefix, "NumEventsDump", num);608 SetDumpEvents(GetEnvValue(env, prefix, "NumEventsDump", fNumEventsDump)); 627 609 rc = kTRUE; 628 610 } 629 SetDumpEvents(num); 630 631 Byte_t max = fMaxHiGainVar; 611 632 612 if (IsEnvDefined(env, prefix, "MaxHiGainVar", print)) 633 613 { 634 max = GetEnvValue(env, prefix, "MaxHiGainVar", max);614 SetMaxHiGainVar(GetEnvValue(env, prefix, "MaxHiGainVar", fMaxHiGainVar)); 635 615 rc = kTRUE; 636 616 } 637 SetMaxHiGainVar(max); 638 639 Bool_t upd = fPedestalUpdate; 617 640 618 if (IsEnvDefined(env, prefix, "PedestalUpdate", print)) 641 619 { 642 upd = GetEnvValue(env, prefix, "PedestalUpdate", upd);620 SetPedestalUpdate(GetEnvValue(env, prefix, "PedestalUpdate", fPedestalUpdate)); 643 621 rc = kTRUE; 644 622 } 645 SetPedestalUpdate(upd);646 623 647 624 return rc; -
trunk/MagicSoft/Mars/mpedestal/MPedCalcFromLoGain.h
r4601 r4609 62 62 63 63 // Setter 64 void SetRange(Byte_t hifirst= 0, Byte_t hilast=0, Byte_t lofirst=0, Byte_t lolast=0);65 void SetWindowSize(Byte_t windowh= 0, Byte_t windowl=0);66 void SetMaxHiGainVar(Byte_t maxvar= 0);64 void SetRange(Byte_t hifirst=fgHiGainFirst, Byte_t hilast=fgHiGainLast, Byte_t lofirst=fgLoGainFirst, Byte_t lolast=fgLoGainLast); 65 void SetWindowSize(Byte_t windowh=fgHiGainWindowSize, Byte_t windowl=fgLoGainWindowSize); 66 void SetMaxHiGainVar(Byte_t maxvar=fgMaxHiGainVar); 67 67 void SetDumpEvents(UInt_t dumpevents = 0) {fNumEventsDump = dumpevents;} 68 68 void SetPedestalUpdate(Bool_t pedupdate) {fPedestalUpdate = pedupdate;} -
trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.cc
r4601 r4609 394 394 395 395 const Int_t npixels = fPedestals->GetSize(); 396 const Int_t areas = fPedestals->GetAverageAreas();397 const Int_t sectors = fPedestals->GetAverageSectors();396 // const Int_t areas = fPedestals->GetAverageAreas(); 397 // const Int_t sectors = fPedestals->GetAverageSectors(); 398 398 399 399 if (fSumx.GetSize()==0) … … 402 402 fSumx2.Set(npixels); 403 403 404 fAreaSumx. Set(areas);405 fAreaSumx2.Set(areas);406 fAreaValid.Set(areas);407 408 fSectorSumx. Set(sectors);409 fSectorSumx2.Set(sectors);410 fSectorValid.Set(sectors);404 // fAreaSumx. Set(areas); 405 // fAreaSumx2.Set(areas); 406 // fAreaValid.Set(areas); 407 408 // fSectorSumx. Set(sectors); 409 // fSectorSumx2.Set(sectors); 410 // fSectorValid.Set(sectors); 411 411 412 412 fSumx.Reset(); … … 439 439 Int_t MPedCalcPedRun::Process() 440 440 { 441 MRawEvtPixelIter pixel(fRawEvt); 442 443 while (pixel.Next()) 444 { 445 const UInt_t idx = pixel.GetPixelId(); 446 447 Byte_t *ptr = pixel.GetHiGainSamples() + fHiGainFirst; 448 Byte_t *end = ptr + fWindowSizeHiGain; 449 450 UInt_t sum = 0; 451 UInt_t sqr = 0; 452 453 if (fWindowSizeHiGain != 0) 454 { 455 do 456 { 457 sum += *ptr; 458 sqr += *ptr * *ptr; 459 } 460 while (++ptr != end); 461 } 462 463 if (fWindowSizeLoGain != 0) 464 { 465 ptr = pixel.GetLoGainSamples() + fLoGainFirst; 466 end = ptr + fWindowSizeLoGain; 467 468 do 469 { 470 sum += *ptr; 471 sqr += *ptr * *ptr; 472 } 473 while (++ptr != end); 474 } 475 476 const Float_t msum = (Float_t)sum; 477 478 // 479 // These three lines have been uncommented by Markus Gaug 480 // If anybody needs them, please contact me!! 481 // 482 // const Float_t higainped = msum/fNumHiGainSlices; 483 // const Float_t higainrms = TMath::Sqrt((msqr-msum*msum/fNumHiGainSlices)/(fNumHiGainSlices-1.)); 484 // (*fPedestals)[idx].Set(higainped, higainrms); 485 486 fSumx[idx] += msum; 487 // 488 // The old version: 489 // 490 // const Float_t msqr = (Float_t)sqr; 491 // fSumx2[idx] += msqr; 492 // 493 // The new version: 494 // 495 const Float_t sqrsum = msum*msum; 496 fSumx2[idx] += sqrsum; 497 } 498 499 fNumSamplesTot += fWindowSizeHiGain + fWindowSizeLoGain; 500 501 return kTRUE; 502 /* 441 503 MRawEvtPixelIter pixel(fRawEvt); 442 504 … … 445 507 const UInt_t idx = pixel.GetPixelId(); 446 508 const UInt_t aidx = (*fGeom)[idx].GetAidx(); 447 const UInt_t sector = (*fGeom)[idx].GetSector(); 509 const UInt_t sector = (*fGeom)[idx].GetSector(); 448 510 449 511 Byte_t *ptr = pixel.GetHiGainSamples() + fHiGainFirst; … … 490 552 fSumx[idx] += msum; 491 553 fAreaSumx[aidx] += msum; 492 fSectorSumx[sector] += msum; 554 fSectorSumx[sector] += msum; 493 555 // 494 556 // The old version: … … 502 564 fSumx2[idx] += sqrsum; 503 565 fAreaSumx2[aidx] += sqrsum; 504 fSectorSumx2[sector] += sqrsum; 566 fSectorSumx2[sector] += sqrsum; 505 567 } 506 568 507 569 fPedestals->SetReadyToSave(); 508 570 fNumSamplesTot += fWindowSizeHiGain + fWindowSizeLoGain; 509 571 */ 510 572 return kTRUE; 511 573 } … … 517 579 Int_t MPedCalcPedRun::PostProcess() 518 580 { 519 // Compute pedestals and rms from the whole run 520 const ULong_t n = fNumSamplesTot; 521 const ULong_t nevts = GetNumExecutions(); 522 523 MRawEvtPixelIter pixel(fRawEvt); 524 525 while (pixel.Next()) 526 { 527 528 const Int_t pixid = pixel.GetPixelId(); 529 const UInt_t aidx = (*fGeom)[pixid].GetAidx(); 530 const UInt_t sector = (*fGeom)[pixid].GetSector(); 531 532 fAreaValid [aidx]++; 533 fSectorValid[sector]++; 534 535 const Float_t sum = fSumx.At(pixid); 536 const Float_t sum2 = fSumx2.At(pixid); 537 const Float_t higainped = sum/n; 538 // 539 // The old version: 540 // 541 // const Float_t higainrms = TMath::Sqrt((sum2-sum*sum/n)/(n-1.)); 542 // 543 // The new version: 544 // 545 // 1. Calculate the Variance of the sums: 546 Float_t higainVar = (sum2-sum*sum/nevts)/(nevts-1.); 547 // 2. Scale the variance to the number of slices: 548 higainVar /= (Float_t)(fWindowSizeHiGain+fWindowSizeLoGain); 549 // 3. Calculate the RMS from the Variance: 550 (*fPedestals)[pixid].Set(higainped, higainVar < 0 ? 0. : TMath::Sqrt(higainVar)); 551 552 } 581 const ULong_t nevts = GetNumExecutions(); 582 if (nevts<1) 583 { 584 *fLog << err << "ERROR - Not enough events recorded...abort." << endl; 585 return kFALSE; 586 } 587 588 *fLog << flush << inf << "Calculating pedestals..." << flush; 589 590 // Compute pedestals and rms from the whole run 591 const ULong_t n = fNumSamplesTot; 592 const ULong_t npix = fGeom->GetNumPixels(); 593 594 for (UInt_t pixidx=0; pixidx<npix; pixidx++) 595 { 596 const Float_t sum = fSumx.At(pixidx); 597 const Float_t sum2 = fSumx2.At(pixidx); 598 const Float_t higainped = sum/n; 599 // 600 // The old version: 601 // 602 // const Float_t higainrms = TMath::Sqrt((sum2-sum*sum/n)/(n-1.)); 603 // 604 // The new version: 605 // 606 // 1. Calculate the Variance of the sums: 607 Float_t higainVar = (sum2-sum*sum/nevts)/(nevts-1.); 608 // 2. Scale the variance to the number of slices: 609 higainVar /= (Float_t)(fWindowSizeHiGain+fWindowSizeLoGain); 610 // 3. Calculate the RMS from the Variance: 611 (*fPedestals)[pixidx].Set(higainped, higainVar < 0 ? 0. : TMath::Sqrt(higainVar), 0, nevts); 612 } 613 614 *fLog << flush << inf << "Calculating means..." << flush; 615 616 fPedestals->SetTotalEntries(fNumSamplesTot); 617 fPedestals->ReCalc(*fGeom); 618 fPedestals->SetReadyToSave(); 619 620 return kTRUE; 553 621 554 622 // 555 623 // Loop over the (two) area indices to get the averaged pedestal per aidx 556 624 // 625 /* 626 627 // Compute pedestals and rms from the whole run 628 const ULong_t n = fNumSamplesTot; 629 const ULong_t nevts = GetNumExecutions(); 630 631 MRawEvtPixelIter pixel(fRawEvt); 632 633 while (pixel.Next()) 634 { 635 const Int_t pixid = pixel.GetPixelId(); 636 const UInt_t aidx = (*fGeom)[pixid].GetAidx(); 637 const UInt_t sector = (*fGeom)[pixid].GetSector(); 638 639 fAreaValid [aidx]++; 640 fSectorValid[sector]++; 641 642 const Float_t sum = fSumx.At(pixid); 643 const Float_t sum2 = fSumx2.At(pixid); 644 const Float_t higainped = sum/n; 645 // 646 // The old version: 647 // 648 // const Float_t higainrms = TMath::Sqrt((sum2-sum*sum/n)/(n-1.)); 649 // 650 // The new version: 651 // 652 // 1. Calculate the Variance of the sums: 653 Float_t higainVar = (sum2-sum*sum/nevts)/(nevts-1.); 654 // 2. Scale the variance to the number of slices: 655 higainVar /= (Float_t)(fWindowSizeHiGain+fWindowSizeLoGain); 656 // 3. Calculate the RMS from the Variance: 657 (*fPedestals)[pixid].Set(higainped, higainVar < 0 ? 0. : TMath::Sqrt(higainVar), 0, nevts); 658 } 557 659 for (Int_t aidx=0; aidx<fAreaValid.GetSize(); aidx++) 558 660 { … … 582 684 fPedestals->GetAverageArea(aidx).Set(higainped, higainrms); 583 685 } 584 686 */ 585 687 // 586 688 // Loop over the (six) sector indices to get the averaged pedestal per sector 587 689 // 690 /* 588 691 for (Int_t sector=0; sector<fSectorValid.GetSize(); sector++) 589 692 { … … 613 716 fPedestals->GetAverageSector(sector).Set(higainped, higainrms); 614 717 } 615 616 718 fPedestals->SetTotalEntries(fNumSamplesTot); 617 719 fPedestals->SetReadyToSave(); 618 720 619 721 return kTRUE; 722 */ 620 723 } 621 724 -
trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.h
r4601 r4609 33 33 TArrayD fSumx; // sum of values 34 34 TArrayD fSumx2; // sum of squared values 35 /* 35 36 TArrayD fAreaSumx; // averaged sum of values per area idx 36 37 TArrayD fAreaSumx2; // averaged sum of squared values per area idx … … 39 40 TArrayD fSectorSumx2; // averaged sum of squared values per sector 40 41 TArrayI fSectorValid; // number of valid pixel with sector idx 42 */ 41 43 42 44 Int_t PreProcess (MParList *pList); -
trunk/MagicSoft/Mars/mpedestal/MPedPhotCam.cc
r4384 r4609 229 229 void MPedPhotCam::ReCalc(const MGeomCam &geom, MBadPixelsCam *bad) 230 230 { 231 232 231 const Int_t np = GetSize(); 233 232 const Int_t ns = GetNumSectors(); 234 233 const Int_t na = GetNumAreas(); 235 236 237 234 238 235 TArrayI acnt(na); … … 243 240 TArrayD ssumr(ns); 244 241 245 246 242 for (int i=0; i<np; i++) 247 243 { 248 249 250 if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun)) 251 continue; //was: .IsBad() 244 if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun)) 245 continue; //was: .IsBad() 252 246 253 247 // Create sums for areas and sectors … … 268 262 ssumr[sect] += ne*rms; 269 263 scnt[sect] += ne; 270 271 272 }273 274 for (int i=0; i<ns; i++){275 if (scnt[i]>0) GetSector(i).Set(ssumx[i]/scnt[i], ssumr[i]/scnt[i], scnt[i]);276 else GetSector(i).Set(-1., -1., 0);277 } 278 279 for (int i=0; i<na; i++){280 if (acnt[i]>0)GetArea(i).Set(asumx[i]/acnt[i], asumr[i]/acnt[i], acnt[i]);281 else GetArea(i).Set(-1., -1., 0);282 }264 } 265 266 for (int i=0; i<ns; i++) 267 if (scnt[i]>0) 268 GetSector(i).Set(ssumx[i]/scnt[i], ssumr[i]/scnt[i], scnt[i]); 269 else 270 GetSector(i).Clear(); 271 272 for (int i=0; i<na; i++) 273 if (acnt[i]>0) 274 GetArea(i).Set(asumx[i]/acnt[i], asumr[i]/acnt[i], acnt[i]); 275 else 276 GetArea(i).Clear(); 283 277 } 284 278 -
trunk/MagicSoft/Mars/mpedestal/MPedPhotCam.h
r4384 r4609 51 51 void Print(Option_t *o="") const; 52 52 53 void ReCalc(const MGeomCam &geom, MBadPixelsCam *bad );53 void ReCalc(const MGeomCam &geom, MBadPixelsCam *bad=NULL); 54 54 55 55 Bool_t GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type=0) const; -
trunk/MagicSoft/Mars/mpedestal/MPedestalCam.cc
r4404 r4609 35 35 #include "MPedestalPix.h" 36 36 37 #include <TArrayI.h> 38 #include <TArrayD.h> 39 37 40 #include <TClonesArray.h> 38 41 … … 43 46 44 47 #include "MGeomCam.h" 48 #include "MGeomPix.h" 49 50 #include "MBadPixelsCam.h" 51 #include "MBadPixelsPix.h" 45 52 46 53 ClassImp(MPedestalCam); … … 311 318 } 312 319 313 Bool_t MPedestalCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const { 314 315 if (GetSize() <= idx) 316 return kFALSE; 317 318 if (!(*this)[idx].IsValid()) 319 return kFALSE; 320 321 switch (type) { 320 // -------------------------------------------------------------------------- 321 // 322 // Calculates the avarage pedestal and pedestal rms for all sectors 323 // and pixel sizes. The geometry container is used to get the necessary 324 // geometry information (sector number, size index) If the bad pixel 325 // container is given all pixels which have the flag 'bad' are ignored 326 // in the calculation of the sector and size average. 327 // 328 void MPedestalCam::ReCalc(const MGeomCam &geom, MBadPixelsCam *bad) 329 { 330 const Int_t np = GetSize(); 331 const Int_t ns = geom.GetNumSectors(); 332 const Int_t na = geom.GetNumAreas(); 333 334 TArrayI acnt(na); 335 TArrayI scnt(ns); 336 TArrayD asumx(na); 337 TArrayD ssumx(ns); 338 TArrayD asumr(na); 339 TArrayD ssumr(ns); 340 341 for (int i=0; i<np; i++) 342 { 343 if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun)) 344 continue; //was: .IsBad() 345 346 // Create sums for areas and sectors 347 const UInt_t aidx = geom[i].GetAidx(); 348 const UInt_t sect = geom[i].GetSector(); 349 350 const MPedestalPix &pix = (*this)[i]; 351 352 const UInt_t ne = pix.GetNumEvents(); 353 const Float_t mean = pix.GetPedestal(); 354 const Float_t rms = pix.GetPedestalRms(); 355 356 asumx[aidx] += ne*mean; 357 asumr[aidx] += ne*rms; 358 acnt[aidx] += ne; 359 360 ssumx[sect] += ne*mean; 361 ssumr[sect] += ne*rms; 362 scnt[sect] += ne; 363 } 364 365 for (int i=0; i<ns; i++) 366 if (scnt[i]>0) 367 GetAverageSector(i).Set(ssumx[i]/scnt[i], ssumr[i]/scnt[i], 0, scnt[i]); 368 else 369 GetAverageSector(i).Clear(); 370 371 for (int i=0; i<na; i++) 372 if (acnt[i]>0) 373 GetAverageArea(i).Set(asumx[i]/acnt[i], asumr[i]/acnt[i], 0, acnt[i]); 374 else 375 GetAverageArea(i).Clear(); 376 } 377 378 Bool_t MPedestalCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const 379 { 380 if (GetSize() <= idx) 381 return kFALSE; 382 383 if (!(*this)[idx].IsValid()) 384 return kFALSE; 385 386 switch (type) 387 { 322 388 case 0: 323 val = (*this)[idx].GetPedestal();324 break;389 val = (*this)[idx].GetPedestal(); 390 break; 325 391 case 1: 326 val = fTotalEntries > 0 ?327 (*this)[idx].GetPedestalRms()/TMath::Sqrt((Float_t)fTotalEntries) 328 329 break;392 val = fTotalEntries > 0 ? 393 (*this)[idx].GetPedestalRms()/TMath::Sqrt((Float_t)fTotalEntries) 394 : (*this)[idx].GetPedestalError(); 395 break; 330 396 case 2: 331 val = (*this)[idx].GetPedestalRms();332 break;397 val = (*this)[idx].GetPedestalRms(); 398 break; 333 399 case 3: 334 val = fTotalEntries > 0 ?335 (*this)[idx].GetPedestalRms()/TMath::Sqrt((Float_t)fTotalEntries)/2. 336 337 break;400 val = fTotalEntries > 0 ? 401 (*this)[idx].GetPedestalRms()/TMath::Sqrt((Float_t)fTotalEntries)/2. 402 : (*this)[idx].GetPedestalRmsError(); 403 break; 338 404 default: 339 return kFALSE;405 return kFALSE; 340 406 } 341 return kTRUE;407 return kTRUE; 342 408 } 343 409 344 410 void MPedestalCam::DrawPixelContent(Int_t idx) const 345 411 { 346 *fLog << warn << "MPedestalCam::DrawPixelContent - not available." << endl;347 } 412 *fLog << warn << "MPedestalCam::DrawPixelContent - not available." << endl; 413 } -
trunk/MagicSoft/Mars/mpedestal/MPedestalCam.h
r3803 r4609 13 13 class MGeomCam; 14 14 class MPedestalPix; 15 class MBadPixelsCam; 15 16 16 17 class MPedestalCam : public MParContainer, public MCamEvent … … 51 52 void InitAverageSectors ( const UInt_t i ); 52 53 54 void ReCalc(const MGeomCam &geom, MBadPixelsCam *bad=NULL); 55 53 56 void Print(Option_t *o="") const; 54 57 -
trunk/MagicSoft/Mars/mpedestal/MPedestalPix.h
r4556 r4609 37 37 Float_t GetPedestalRmsError() const { return fNumEvents>0 ? fPedestalRms/TMath::Sqrt((Float_t)fNumEvents/2) : 0; } 38 38 39 UInt_t GetNumEvents() const { return fNumEvents; } 40 39 41 Bool_t IsValid() const; 40 42 -
trunk/MagicSoft/Mars/mraw/MRawEvtHeader.cc
r4577 r4609 74 74 // 75 75 // UInt_t fNumTrigLvl2 76 // ------------------ 76 // -------------------- 77 77 // 78 78 // Number of second level trigger … … 90 90 // Type of Trigger. 91 91 // This is a Byte (8 bit) to indicated if any of the pixels 92 // have a non-negligible low gain (1) or not (0) 92 // have a non-negligible low gain (1) or not (0) 93 // 94 // UInt_t fCalibPattern 95 // -------------------- 96 // 97 // Bits 1-16: Pulser slot pattern: 16 LEDs slots (see Diploma Thesis 98 // of Michele) 99 // 100 // Bits 17-24: Power of Continous light source: 256 level 101 // 102 // Bits 25-28: Farbe der Continuous Light: red-amber-green/blue-uv 103 // 104 // Bit 29: Calibration trigger On/off 105 // Bit 30: Pedestal trigger on/off 106 // Bit 31: PIN Diode calibration trigger on/off 107 // 108 // 109 // Class Version 2: 110 // --------------- 111 // - added fCalibPattern 93 112 // 94 113 ///////////////////////////////////////////////////////////////////////////// … … 379 398 return fTrigPattern[0]; 380 399 } 400 401 // -------------------------------------------------------------------------- 402 // 403 // return pulser slot patter, see class reference 404 // 405 UShort_t MRawEvtHeader::GetPulserSlotPattern() const 406 { 407 return fTrigPattern[1] & 0xffff; 408 } 409 410 // -------------------------------------------------------------------------- 411 // 412 // return power of continous light, see class reference 413 // 414 Byte_t MRawEvtHeader::GetPowerOfContLight() const 415 { 416 return (fTrigPattern[1]<<16) & 0xff; 417 } 418 419 // -------------------------------------------------------------------------- 420 // 421 // return color of continous light, see class reference 422 // 423 MRawEvtHeader::CalibCol_t MRawEvtHeader::GetContLedColor() const 424 { 425 return (CalibCol_t)((fTrigPattern[1]<<24)&0xf); 426 } -
trunk/MagicSoft/Mars/mraw/MRawEvtHeader.h
r4577 r4609 24 24 kTTPedestal = 1, 25 25 kTTCalibration = 2 26 }; 27 28 enum CalibCol_t { 29 kColRed = BIT(0), 30 kColAmber = BIT(1), 31 kColGreen = BIT(2), 32 kColBlue = BIT(3), 33 kColUV = BIT(4) 26 34 }; 27 35 … … 60 68 void FillHeader(UInt_t, Float_t=0); 61 69 62 UShort_t GetTrigType() const { return fTrigType; }63 UInt_t GetNumTrigLvl1() const { return fNumTrigLvl1; }64 UInt_t GetNumTrigLvl2() const { return fNumTrigLvl2; }65 UInt_t GetDAQEvtNumber() const { return fDAQEvtNumber; }70 UShort_t GetTrigType() const { return fTrigType; } 71 UInt_t GetNumTrigLvl1() const { return fNumTrigLvl1; } 72 UInt_t GetNumTrigLvl2() const { return fNumTrigLvl2; } 73 UInt_t GetDAQEvtNumber() const { return fDAQEvtNumber; } 66 74 67 UInt_t GetTriggerID() const; 75 UInt_t GetTriggerID() const; 76 77 UShort_t GetPulserSlotPattern() const; 78 Byte_t GetPowerOfContLight() const; 79 CalibCol_t GetContLedColor() const; 68 80 69 81 Int_t ReadEvt(istream& fin, UShort_t ver); 70 void SkipEvt(istream& fin, UShort_t ver);82 void SkipEvt(istream& fin, UShort_t ver); 71 83 72 84 ClassDef(MRawEvtHeader, 1) // Parameter Conatiner for raw EVENT HEADER
Note:
See TracChangeset
for help on using the changeset viewer.