Changeset 4900
- Timestamp:
- 09/09/04 17:23:23 (20 years ago)
- Location:
- trunk/MagicSoft/Mars
- Files:
-
- 1 added
- 5 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/MagicSoft/Mars/Changelog
r4899 r4900 38 38 and leave only typical Gaussian behaviour members in MHGausEvents. 39 39 - all MHCalibration*Pix classes derive now from MHCalibrationPix 40 41 * mcalib/MHCalibrationCam.[h,cc] 42 - operators return MHCalibrationPix instead of MHGausEvents 43 - new function ResetHists() needed by MCalibColorSteer 44 - removed TArrayI *Overflow. 40 45 41 46 * mcalib/MCalibColorSet.[h,cc] -
trunk/MagicSoft/Mars/mcalib/CalibLinkDef.h
r4858 r4900 6 6 7 7 #pragma link C++ class MCalibColorSet+; 8 #pragma link C++ class MCalibColorSteer+; 8 9 9 10 #pragma link C++ class MCalibrate+; … … 11 12 #pragma link C++ class MCalibrateRelTimes+; 12 13 14 #pragma link C++ class MCalibrationIntensityCam+; 15 #pragma link C++ class MCalibrationIntensityChargeCam+; 16 #pragma link C++ class MCalibrationIntensityQECam+; 17 #pragma link C++ class MCalibrationIntensityRelTimeCam+; 13 18 #pragma link C++ class MCalibrationCam+; 14 19 #pragma link C++ class MCalibrationPix+; … … 33 38 34 39 #pragma link C++ class MHCalibrationCam+; 40 #pragma link C++ class MHCalibrationPix+; 35 41 #pragma link C++ class MHCalibrationChargeCam+; 36 42 #pragma link C++ class MHCalibrationChargePix+; -
trunk/MagicSoft/Mars/mcalib/MHGausEvents.cc
r4873 r4900 101 101 using namespace std; 102 102 103 const Int_t MHGausEvents::fgBinsAfterStripping = 40;104 const Float_t MHGausEvents::fgBlackoutLimit = 5.;105 103 const Int_t MHGausEvents::fgNDFLimit = 2; 106 const Float_t MHGausEvents::fgPickupLimit = 5.;107 104 const Float_t MHGausEvents::fgProbLimit = 0.001; 108 const Int_t MHGausEvents::fgPowerProbabilityBins = 20;105 const Int_t MHGausEvents::fgPowerProbabilityBins = 30; 109 106 // -------------------------------------------------------------------------- 110 107 // … … 114 111 // - the default Probability Limit for fProbLimit (fgProbLimit) 115 112 // - the default NDF Limit for fNDFLimit (fgNDFLimit) 116 // - the default number for fPickupLimit (fgPickupLimit)117 // - the default number for fBlackoutLimit (fgBlackoutLimit)118 113 // - the default number of bins after stripping for fBinsAfterStipping (fgBinsAfterStipping) 119 114 // - the default name of the fHGausHist ("HGausHist") … … 132 127 // 133 128 MHGausEvents::MHGausEvents(const char *name, const char *title) 134 : f EventFrequency(0), fHPowerProbability(NULL),129 : fHPowerProbability(NULL), 135 130 fPowerSpectrum(NULL), 131 fFGausFit(NULL), fFExpFit(NULL), 132 fFirst(0.), 136 133 fGraphEvents(NULL), fGraphPowerSpectrum(NULL), 137 fFGausFit(NULL), fFExpFit(NULL), 138 fFirst(0.), fLast(100.), 139 fNbins(100), fPixId(-1) 134 fLast(100.), fNbins(100) 140 135 { 141 136 … … 146 141 147 142 SetBinsAfterStripping(); 148 SetBlackoutLimit();149 143 SetNDFLimit(); 150 SetPickupLimit();151 144 SetPowerProbabilityBins(); 152 145 SetProbLimit(); … … 206 199 // Sets: 207 200 // - all other pointers to NULL 208 // - all variables to 0. , except fPixId to -1 and keep fEventFrequency201 // - all variables to 0. 209 202 // - all flags to kFALSE 210 203 // … … 226 219 fProb = 0.; 227 220 228 fSaturated = 0;229 221 fCurrentSize = 0; 230 222 … … 271 263 272 264 273 // -----------------------------------------------------------------------------274 //275 // Bypasses the Gauss fit by taking mean and RMS from the histogram276 //277 // Errors are determined in the following way:278 // MeanErr = RMS / Sqrt(entries)279 // SigmaErr = RMS / (2.*Sqrt(entries) )280 //281 void MHGausEvents::BypassFit()282 {283 284 const Stat_t entries = fHGausHist.GetEntries();285 286 if (entries <= 0.)287 {288 *fLog << warn << GetDescriptor()289 << ": Cannot bypass fit. Number of entries smaller or equal 0 in pixel: " << fPixId << endl;290 return;291 }292 293 fMean = fHGausHist.GetMean();294 fMeanErr = fHGausHist.GetRMS() / TMath::Sqrt(entries);295 fSigma = fHGausHist.GetRMS() ;296 fSigmaErr = fHGausHist.GetRMS() / TMath::Sqrt(entries) / 2.;297 }298 299 300 301 // --------------------------------------------------------------------------302 //303 // - Set fPixId to id304 //305 // Add id to names and titles of:306 // - fHGausHist307 //308 void MHGausEvents::ChangeHistId(const Int_t id)309 {310 311 fPixId = id;312 313 fHGausHist.SetName( Form("%s%d", fHGausHist.GetName(), id));314 fHGausHist.SetTitle( Form("%s%d", fHGausHist.GetTitle(), id));315 316 fName = Form("%s%d", fName.Data(), id);317 fTitle = Form("%s%d", fTitle.Data(), id);318 319 }320 265 321 266 // ----------------------------------------------------------------------------- … … 328 273 Float_t *xaxis = new Float_t[n]; 329 274 330 if (fEventFrequency) 331 for (Int_t i=0;i<n;i++) 332 xaxis[i] = (Float_t)i/fEventFrequency; 333 else 334 for (Int_t i=0;i<n;i++) 335 xaxis[i] = (Float_t)i; 336 275 for (Int_t i=0;i<n;i++) 276 xaxis[i] = (Float_t)i; 277 337 278 return xaxis; 338 279 … … 344 285 // Create the fourier spectrum using the class MFFT. 345 286 // The result is projected into a histogram and fitted by an exponential 346 // fHPowerProbability->SetDirectory(NULL);287 // 347 288 void MHGausEvents::CreateFourierSpectrum() 348 289 { … … 353 294 if (fEvents.GetSize() < 8) 354 295 { 355 *fLog << warn << "Cannot create Fourier spectrum in pixel: " << fPixId296 *fLog << warn << "Cannot create Fourier spectrum in: " << fName 356 297 << ". Number of events smaller than 8 " << endl; 357 298 return; 358 299 } 359 300 360 361 301 // 362 302 // The number of entries HAS to be a potence of 2, … … 436 376 const Int_t n = fEvents.GetSize(); 437 377 378 if (n==0) 379 return; 380 438 381 fGraphEvents = new TGraph(n,CreateEventXaxis(n),fEvents.GetArray()); 439 382 fGraphEvents->SetTitle("Evolution of Events with time"); 440 fGraphEvents->GetXaxis()->SetTitle((fEventFrequency) ? "Time [s]" : "Event Nr.");441 383 fGraphEvents->GetYaxis()->SetTitle(fHGausHist.GetXaxis()->GetTitle()); 442 384 fGraphEvents->GetYaxis()->CenterTitle(); … … 457 399 fGraphPowerSpectrum = new TGraph(n,CreatePSDXaxis(n),fPowerSpectrum->GetArray()); 458 400 fGraphPowerSpectrum->SetTitle("Power Spectrum Density"); 459 fGraphPowerSpectrum->GetXaxis()->SetTitle((fEventFrequency) ? "Frequency [Hz]" : "Frequency");460 401 fGraphPowerSpectrum->GetYaxis()->SetTitle("P(f)"); 461 402 fGraphPowerSpectrum->GetYaxis()->CenterTitle(); … … 473 414 Float_t *xaxis = new Float_t[n]; 474 415 475 if (fEventFrequency) 476 for (Int_t i=0;i<n;i++) 477 xaxis[i] = 0.5*(Float_t)i*fEventFrequency/n; 478 else 479 for (Int_t i=0;i<n;i++) 480 xaxis[i] = 0.5*(Float_t)i/n; 481 416 for (Int_t i=0;i<n;i++) 417 xaxis[i] = 0.5*(Float_t)i/n; 418 482 419 return xaxis; 483 420 … … 534 471 pad->cd(1); 535 472 536 if (!IsEmpty() )473 if (!IsEmpty() && !IsOnlyOverflow() && !IsOnlyUnderflow()) 537 474 gPad->SetLogy(); 538 475 … … 571 508 CreateGraphEvents(); 572 509 510 if (!fGraphEvents) 511 return; 512 573 513 fGraphEvents->SetBit(kCanDelete); 574 514 fGraphEvents->SetTitle("Events with time"); … … 666 606 // 667 607 // First, strip the zeros from the edges which contain only zeros and rebin 668 // to about fBinsAfterStripping bins. 608 // to fBinsAfterStripping bins. If fBinsAfterStripping is 0, reduce bins only by 609 // a factor 10. 669 610 // 670 611 // (ATTENTION: The Chisquare method is more sensitive, 671 612 // the _less_ bins, you have!) 672 613 // 673 StripZeros(&fHGausHist,fBinsAfterStripping); 614 StripZeros(&fHGausHist, 615 fBinsAfterStripping ? fBinsAfterStripping 616 : (fNbins > 1000 ? fNbins/10 : 0)); 674 617 675 618 TAxis *axe = fHGausHist.GetXaxis(); … … 693 636 { 694 637 *fLog << warn << dbginf << "WARNING: Could not create fit function for Gauss fit " 695 << "in pixel: " << fPixId<< endl;638 << "in: " << fName << endl; 696 639 return kFALSE; 697 640 } … … 731 674 } 732 675 733 // -------------------------------------------------------------------------------734 //735 // Return the number of "blackout" events, which are events with values higher736 // than fBlackoutLimit sigmas from the mean737 //738 //739 const Double_t MHGausEvents::GetBlackout() const740 {741 742 if ((fMean == 0.) && (fSigma == 0.))743 return -1.;744 745 const Int_t first = fHGausHist.GetXaxis()->GetFirst();746 const Int_t last = fHGausHist.GetXaxis()->FindBin(fMean-fBlackoutLimit*fSigma);747 748 if (first >= last)749 return 0.;750 751 return fHGausHist.Integral(first, last, "width");752 }753 676 754 677 const Double_t MHGausEvents::GetChiSquare() const … … 786 709 } 787 710 788 789 // -------------------------------------------------------------------------------790 //791 // Return the number of "pickup" events, which are events with values higher792 // than fPickupLimit sigmas from the mean793 //794 //795 const Double_t MHGausEvents::GetPickup() const796 {797 798 if ((fMean == 0.) && (fSigma == 0.))799 return -1.;800 801 const Int_t first = fHGausHist.GetXaxis()->FindBin(fMean+fPickupLimit*fSigma);802 const Int_t last = fHGausHist.GetXaxis()->GetLast();803 804 if (first >= last)805 return 0.;806 807 return fHGausHist.Integral(first, last, "width");808 }809 711 810 712 … … 890 792 891 793 *fLog << all << endl; 892 *fLog << all << "Results of the Gauss Fit in pixel: " << fPixId<< endl;794 *fLog << all << "Results of the Gauss Fit in: " << fName << endl; 893 795 *fLog << all << "Mean: " << GetMean() << endl; 894 796 *fLog << all << "Sigma: " << GetSigma() << endl; … … 900 802 } 901 803 902 // --------------------------------------------------------------------------903 //904 // Re-normalize the results, has to be overloaded905 //906 void MHGausEvents::Renorm()907 {908 }909 910 // -----------------------------------------------------------------------------911 //912 // If flag IsGausFitOK() is set (histogram already successfully fitted),913 // returns kTRUE914 //915 // If both fMean and fSigma are still zero, call FitGaus()916 //917 // Repeats the Gauss fit in a smaller range, defined by:918 //919 // min = GetMean() - fBlackoutLimit * GetSigma();920 // max = GetMean() + fPickupLimit * GetSigma();921 //922 // The fit results are retrieved and stored in class-own variables.923 //924 // A flag IsGausFitOK() is set according to whether the fit probability925 // is smaller or bigger than fProbLimit, whether the NDF is bigger than926 // fNDFLimit and whether results are NaNs.927 //928 Bool_t MHGausEvents::RepeatFit(const Option_t *option)929 {930 931 if (IsGausFitOK())932 return kTRUE;933 934 if ((fMean == 0.) && (fSigma == 0.))935 return FitGaus();936 937 //938 // Get new fitting ranges939 //940 Axis_t rmin = fMean - fBlackoutLimit * fSigma;941 Axis_t rmax = fMean + fPickupLimit * fSigma;942 943 Axis_t hmin = fHGausHist.GetBinCenter(fHGausHist.GetXaxis()->GetFirst());944 Axis_t hmax = fHGausHist.GetBinCenter(fHGausHist.GetXaxis()->GetLast()) ;945 946 fFGausFit->SetRange(hmin < rmin ? rmin : hmin , hmax > rmax ? rmax : hmax);947 948 fHGausHist.Fit(fFGausFit,option);949 950 fMean = fFGausFit->GetParameter(1);951 fSigma = fFGausFit->GetParameter(2);952 fMeanErr = fFGausFit->GetParError(1) ;953 fSigmaErr = fFGausFit->GetParError(2) ;954 fProb = fFGausFit->GetProb() ;955 956 //957 // The fit result is accepted under condition:958 // 1) The results are not nan's959 // 2) The NDF is not smaller than fNDFLimit (default: fgNDFLimit)960 // 3) The Probability is greater than fProbLimit (default: fgProbLimit)961 //962 if ( TMath::IsNaN ( fMean )963 || TMath::IsNaN ( fMeanErr )964 || TMath::IsNaN ( fProb )965 || TMath::IsNaN ( fSigma )966 || TMath::IsNaN ( fSigmaErr )967 || fFGausFit->GetNDF() < fNDFLimit968 || fProb < fProbLimit )969 return kFALSE;970 971 SetGausFitOK(kTRUE);972 return kTRUE;973 974 }975 804 976 805 // -------------------------------------------------------------------------- -
trunk/MagicSoft/Mars/mcalib/MHGausEvents.h
r4882 r4900 21 21 private: 22 22 23 const static Int_t fgBinsAfterStripping; //! Default for fBinsAfterStripping (now set to: 40)24 const static Float_t fgBlackoutLimit; //! Default for fBlackoutLimit (now set to: 5. )25 23 const static Int_t fgNDFLimit; //! Default for fNDFLimit (now set to: 2) 26 24 const static Float_t fgProbLimit; //! Default for fProbLimit (now set to: 0.001) 27 const static Float_t fgPickupLimit; //! Default for fPickupLimit (now set to: 5. )28 25 const static Int_t fgPowerProbabilityBins; //! Default for fPowerProbabilityBins (now set to: 20) 29 26 30 27 Int_t fBinsAfterStripping; // Bins for the Gauss Histogram after stripping off the zeros at both ends 31 28 Int_t fCurrentSize; // Current size of the array fEvents 32 Float_t fEventFrequency; // Event frequency in Hertz (to be set)33 29 Byte_t fFlags; // Bit field for the fit result bits 34 30 Int_t fPowerProbabilityBins; // Bins for the projected power spectrum … … 36 32 TH1I *fHPowerProbability; // Fourier transform of fEvents projected on y-axis 37 33 TArrayF *fPowerSpectrum; // Fourier transform of fEvents 38 TGraph *fGraphEvents; //! TGraph to display the event array (will not be cloned!!)39 TGraph *fGraphPowerSpectrum; //! TGraph to display the power spectrum array (will not be cloned!!)40 34 41 enum { kGausFitOK, kExpFitOK, kFourierSpectrumOK, kExcluded }; // Bits for information about fit results 35 enum { kGausFitOK, 36 kExpFitOK, 37 kFourierSpectrumOK, 38 kExcluded }; // Bits for information about fit results 42 39 43 40 protected: 44 41 45 Float_t fBlackoutLimit; // Lower number sigmas from mean until event is considered blackout46 42 TArrayF fEvents; // Array which holds the entries of GausHist 47 43 TF1 *fFGausFit; // Gauss fit for fHGausHist 48 44 TF1 *fFExpFit; // Exponential fit for FHPowerProbability 49 45 Axis_t fFirst; // Lower histogram edge for fHGausHist (used by InitBins()) 46 TGraph *fGraphEvents; //! TGraph to display the event array (will not be cloned!!) 47 TGraph *fGraphPowerSpectrum; //! TGraph to display the power spectrum array (will not be cloned!!) 50 48 TH1F fHGausHist; // Histogram to hold the Gaussian distribution 51 49 Axis_t fLast; // Upper histogram edge for fHGausHist (used by InitBins()) … … 54 52 Int_t fNbins; // Number histogram bins for fHGausHist (used by InitBins()) 55 53 Int_t fNDFLimit; // NDF limit for judgement if fit is OK 56 Int_t fSaturated; // Number of events classified as saturated57 54 Double_t fSigma; // Sigma of the Gauss fit 58 55 Double_t fSigmaErr; // Error of the sigma of the Gauss fit 59 Float_t fPickupLimit; // Upper number sigmas from mean until event is considered pickup60 Int_t fPixId; // Pixel ID61 56 Double_t fProb; // Probability of the Gauss fit 62 57 Float_t fProbLimit; // Probability limit for judgement if fit is OK 63 58 64 Float_t *CreateEventXaxis(Int_t n); // Create an x-axis for the Event TGraphs 65 Float_t *CreatePSDXaxis(Int_t n); // Create an x-axis for the PSD TGraphs 59 virtual Float_t *CreateEventXaxis(Int_t n); // Create an x-axis for the Event TGraphs 60 virtual Float_t *CreatePSDXaxis(Int_t n); // Create an x-axis for the PSD TGraphs 61 virtual void CreateGraphEvents(); // Create the TGraph fGraphEvents of fEvents 62 virtual void CreateGraphPowerSpectrum(); // Create the TGraph fGraphPowerSpectrum out of fPowerSpectrum 66 63 67 64 void DrawPowerSpectrum(TVirtualPad &pad, Int_t i); // Draw graph of fPowerSpectrum and fHPowerProbability 68 65 69 66 // Setters 70 void SetBinsAfterStripping ( const Int_t nbins= fgBinsAfterStripping ){ fBinsAfterStripping =nbins; }67 void SetBinsAfterStripping ( const Int_t nbins=0 ) { fBinsAfterStripping =nbins; } 71 68 void SetPowerProbabilityBins ( const Int_t nbins=fgPowerProbabilityBins ) { fPowerProbabilityBins=nbins; } 72 69 73 70 public: 74 71 75 72 MHGausEvents(const char* name=NULL, const char* title=NULL); … … 85 82 86 83 // Getters 87 const Double_t GetBlackout() const;88 84 const Double_t GetChiSquare() const; 89 85 const Double_t GetExpChiSquare() const; … … 110 106 const Int_t GetNdf() const; 111 107 const Double_t GetOffset() const; 112 const Double_t GetPickup() const;113 const Int_t GetPixId() const { return fPixId; }114 108 TArrayF *GetPowerSpectrum() { return fPowerSpectrum; } 115 109 const TArrayF *GetPowerSpectrum() const { return fPowerSpectrum; } 116 110 const Double_t GetProb() const { return fProb; } 117 const Float_t GetSaturated() const { return fSaturated; }118 111 const Double_t GetSigma() const { return fSigma; } 119 112 const Double_t GetSigmaErr() const { return fSigmaErr; } … … 137 130 const Double_t xmin=0., 138 131 const Double_t xmax=0.); // Fit the histogram HGausHist with a Gaussian 139 Bool_t RepeatFit(const Option_t *option="RQ0"); // Repeat fit within limits defined by fPickupLimit140 void BypassFit(); // Take mean and RMS from the histogram141 132 142 133 // Prints … … 144 135 145 136 // Setters 146 void SetBlackoutLimit ( const Float_t lim=fgBlackoutLimit ) { fBlackoutLimit = lim; }147 void SetEventFrequency ( const Float_t f ) { fEventFrequency = f; }148 137 void SetExcluded ( const Bool_t b=kTRUE ); 149 138 void SetExpFitOK ( const Bool_t b=kTRUE ); … … 156 145 void SetNbins ( const Int_t i ) { fNbins = i; } 157 146 void SetNDFLimit ( const Int_t lim=fgNDFLimit ) { fNDFLimit = lim; } 158 void SetPickupLimit ( const Float_t lim=fgPickupLimit ) { fPickupLimit = lim; }159 void SetPixId ( const Int_t i ) { fPixId = i; }160 147 void SetProb ( const Double_t d ) { fProb = d; } 161 148 void SetProbLimit ( const Float_t lim=fgProbLimit ) { fProbLimit = lim; } 162 void SetSaturated ( const Int_t i ) { fSaturated += i; }163 149 void SetSigma ( const Double_t d ) { fSigma = d; } 164 150 void SetSigmaErr ( const Double_t d ) { fSigmaErr = d; } 165 151 166 // Miscelleaneous167 virtual void ChangeHistId(const Int_t id); // Changes names and titles of the histogram168 virtual void Renorm(); // Re-normalize the results169 170 152 void CreateFourierSpectrum(); // Create the fourier spectrum out of fEvents 171 void CreateGraphEvents(); // Create the TGraph fGraphEvents of fEvents172 void CreateGraphPowerSpectrum(); // Create the TGraph fGraphPowerSpectrum out of fPowerSpectrum173 153 174 154 ClassDef(MHGausEvents, 1) // Base class for events with Gaussian distributed values -
trunk/MagicSoft/Mars/mcalib/Makefile
r4858 r4900 33 33 SRCFILES = MCalibrate.cc \ 34 34 MCalibColorSet.cc \ 35 MCalibColorSteer.cc \ 35 36 MCalibrateData.cc \ 36 37 MCalibrateRelTimes.cc \ 38 MCalibrationIntensityCam.cc \ 39 MCalibrationIntensityChargeCam.cc \ 40 MCalibrationIntensityQECam.cc \ 41 MCalibrationIntensityRelTimeCam.cc \ 37 42 MCalibrationCam.cc \ 38 43 MCalibrationPix.cc \ … … 59 64 MHCalibrationChargePix.cc \ 60 65 MHCalibrationCam.cc \ 66 MHCalibrationPix.cc \ 61 67 MHCalibrationChargeCam.cc \ 62 68 MHCalibrationChargeHiGainPix.cc \
Note:
See TracChangeset
for help on using the changeset viewer.