Index: trunk/MagicSoft/Mars/mtemp/Changelog =================================================================== --- trunk/MagicSoft/Mars/mtemp/Changelog (revision 6446) +++ trunk/MagicSoft/Mars/mtemp/Changelog (revision 6447) @@ -18,4 +18,9 @@ -*-*- END OF LINE -*-*- + +2005/02/14: Hendrik Bartko + * mtemp/mmpi/macros/noise_autocorrelation_AB_all.C + - add new macro to calculate the noise autocorrelation for all + inner pixels @@ -24,5 +29,4 @@ - A wrapper class which applies the general unfolding algorithms of MUnfold to a gamma spectrum - Index: trunk/MagicSoft/Mars/mtemp/mmpi/macros/noise_autocorrelation_AB_all.C =================================================================== --- trunk/MagicSoft/Mars/mtemp/mmpi/macros/noise_autocorrelation_AB_all.C (revision 6447) +++ trunk/MagicSoft/Mars/mtemp/mmpi/macros/noise_autocorrelation_AB_all.C (revision 6447) @@ -0,0 +1,358 @@ +/* ======================================================================== *\ +! +! * +! * This file is part of MARS, the MAGIC Analysis and Reconstruction +! * Software. It is distributed to you in the hope that it can be a useful +! * and timesaving tool in analysing Data of imaging Cerenkov telescopes. +! * It is distributed WITHOUT ANY WARRANTY. +! * +! * Permission to use, copy, modify and distribute this software and its +! * documentation for any purpose is hereby granted without fee, +! * provided that the above copyright notice appear in all copies and +! * that both that copyright notice and this permission notice appear +! * in supporting documentation. It is provided "as is" without express +! * or implied warranty. +! * +! +! +! Author(s): Hendrik Bartko 09/2004 +! +! Copyright: MAGIC Software Development, 2000-2001 +! +! +\* ======================================================================== */ + + +const TString defname = "/.magic/magicserv01/scratch/hbartko/rootdata/2004_09_03/20040903_35698_P_multiplexFADC_E.root"; +const TString defpedname = "/.magic/magicserv01/scratch/hbartko/rootdata/2004_09_03/20040903_35698_P_multiplexFADC_E.root"; + +// works with one pedestal file to calculate the noise autocorrelation + + + + +void noise_autocorrelation_AB_all(const TString fname = defpedname, const TString pedname = defpedname) { + + + // for (int i=1; i++; i<397){ + + TH2F bla; + TH2F * hcorr = new TH2F("hcorr","hcorr",15,0,15,15,0,15); + + + for (int i=1; i<397; i++){ + bla = noise_autocorrelation_AB_pix(fname,pedname, i); + + hcorr->Add(&bla); + } + + + TCanvas * c_corr = new TCanvas("c_corr","c_corr",600,400); + c_corr->SetFillColor(0); + c_corr->SetBorderMode(0); + c_corr->SetGridx(); + c_corr->SetGridy(); + gStyle->SetPalette(1); + + hcorr->Scale(1./396.); + + hcorr->Draw("colz"); + // } + +} + + + + + +TH2F noise_autocorrelation_AB_pix(const TString fname = defpedname, const TString pedname = defpedname, Int_t ipix = 67) { + + MParList plist_ped; + + MTaskList tlist_ped; + plist_ped.AddToList(&tlist_ped); + + MPedestalCam pedcam; + plist_ped.AddToList(&pedcam); + + MReadMarsFile read("Events", pedname); + read.DisableAutoScheme(); + tlist_ped.AddToList(&read); + + MGeomApply geomapl_ped; + MGeomCamMagic geomcam; + tlist_ped.AddToList(&geomapl_ped); + + + MPedCalcFromLoGain pedcalc_ped; + //pedcalc_ped.SetMaxHiGainVar(20); + //pedcalc_ped.SetRange(0, 11, 1, 14); + //pedcalc_ped.SetWindowSize(12,14); + pedcalc_ped.SetPedestalUpdate(kFALSE); + tlist_ped.AddToList(&pedcalc_ped); + + MEvtLoop evtloop_ped; + evtloop_ped.SetParList(&plist_ped); + + if (!evtloop_ped.Eventloop()) + return; + + tlist_ped.PrintStatistics(); + + + // now the event loop for the signal reconstruction with pedestals subtracted + + + MParList plist; + MTaskList tlist; + plist.AddToList(&tlist); + + plist.AddToList(&pedcam); + + MRawRunHeader runheader; + plist.AddToList(&runheader); + + MRawEvtData evtdata; + plist.AddToList(&evtdata); + + MExtractedSignalCam sigcam; + plist.AddToList(&sigcam); + + MArrivalTimeCam timecam; + plist.AddToList(&timecam); + + + MReadMarsFile read("Events", fname); + read.DisableAutoScheme(); + + MGeomApply geomapl; + + Int_t WindowSize = 6; // size of the chosen integration window for signal and time extraction + + MExtractFixedWindowPeakSearch extractor; + extractor.SetRange(0, 14, 0, 14); + extractor.SetWindows( WindowSize, WindowSize, 4); + + MExtractTimeHighestIntegral timeext; + timeext.SetRange(0, 14, 0, 14); + timeext.SetWindowSize(WindowSize,WindowSize); + + + tlist.AddToList(&read); + tlist.AddToList(&geomapl); + tlist.AddToList(&extractor); + tlist.AddToList(&timeext); + + + MEvtLoop evtloop; + evtloop.SetParList(&plist); + + if (!evtloop.PreProcess()) + return; + + TString title = "Noise Autocorrelation B_{ij} = - *, Pixel "; + title += ipix; + + + Int_t First = 1; + + TH2F hcorr("hcorr","hcorr",15,0,15,15,0,15); + + + Float_t B[15][15]; + Float_t bi[15]; + Float_t bj[15]; + + for (int i=0; i<15; i++){ + bi[i]=0; + bj[i]=0; + for (int j=0; j<15; j++){ + B[i][j]=0; + } + } + + const Float_t ped_mean = pedcam[ipix].GetPedestal(); + const Float_t pedRMS = pedcam[ipix].GetPedestalRms(); + + cout << " ped_mean " << ped_mean << " pedRMS " << pedRMS << " ABoffs " << pedcam[ipix].GetPedestalABoffset() << endl; + int count = 0; + + while (tlist.Process()) { + + if (First) { + First = 0; + + const Int_t nh = runheader.GetNumSamplesHiGain(); + const Int_t nl = runheader.GetNumSamplesLoGain(); + const Int_t nt = nh+nl; + + } + + + + MRawEvtPixelIter pixel(&evtdata); + pixel.Jump(ipix); + + Bool_t ABFlag = pixel.HasABFlag(); + + + count++; // couter for the number of events + + const Byte_t *higains = pixel.GetHiGainSamples(); + const Byte_t *logains = pixel.GetLoGainSamples(); + + const Float_t ABoffs = pedcam[ipix].GetPedestalABoffset(); + + Float_t PedMean[2]; + PedMean[0] = ped_mean + ABoffs; + PedMean[1] = ped_mean - ABoffs; + + + // calculate the noise autocorrelation matrix + for (int slice=0; sliceFill(i,j,B[i][j]/count-bi[i]/count*bj[j]/count); + // cout << "i " << i << " j " << j << " B[i][j] " << B[i][j] << " bi[i] " << bi[i] << " bj[j] " << bj[j] << endl; + } + } + + cout << " pedmean " << ped_mean << " pedRMS " << pedRMS << endl; + + + return hcorr; + + /* + + TCanvas * c_corr = new TCanvas("c_corr","c_corr",600,400); + c_corr->SetFillColor(0); + c_corr->SetBorderMode(0); + c_corr->SetGridx(); + c_corr->SetGridy(); + gStyle->SetPalette(1); + + hcorr->SetStats(0); + hcorr->SetTitle(title); + hcorr->SetXTitle("slice i"); + hcorr->SetYTitle("slice j"); + hcorr->Draw("colz"); + + + TCanvas * c_corr_px = new TCanvas("c_corr_px","c_corr_px",600,400); + c_corr_px->SetFillColor(0); + c_corr_px->SetBorderMode(0); + c_corr_px->SetGridx(); + c_corr_px->SetGridy(); + + TH1F * hcorr_px = (TH1F*)hcorr->ProjectionX("hcorr_px",6,6); + hcorr_px->SetLineColor(4); + hcorr_px->SetLineWidth(2); + hcorr_px->SetStats(0); + TString title1 = "Noise Autocorrelation, FADC sample 6, Pixel"; + title1+=ipix; + hcorr_px->SetTitle(title); + hcorr_px->SetXTitle("slice i"); + hcorr_px->SetYTitle("B_{ij} = - * [FADC counts]"); + + hcorr_px->Draw(); + + + + TCanvas * c_shape = new TCanvas("c_shape","c_shape",600,400); + c_shape->SetFillColor(0); + c_shape->SetBorderMode(0); + c_shape->SetGridx(); + c_shape->SetGridy(); + shape->SetStats(0); + TString title2 = "Noise Signal Shape, Signal #geq 3 FADC counts, Pixel "; + title2+=ipix; + shape->SetTitle(title2); + shape->SetXTitle("FADC sample no."); + shape->SetYTitle("signal [FADC counts]"); + shape->SetMarkerStyle(20); + shape->SetMarkerColor(4); + shape->SetMarkerSize(0.3); + shape->Draw(); + + + TCanvas * c_shape_corr = new TCanvas("c_shape_corr","c_shape_corr",600,400); + c_shape_corr->SetFillColor(0); + c_shape_corr->SetBorderMode(0); + c_shape_corr->SetGridx(); + c_shape_corr->SetGridy(); + shape_corr->SetStats(0); + TString title3 = "Noise Signal Shape, Signal #geq 3 FADC counts, Pixel "; + title3+=ipix; + shape_corr->SetTitle(title3); + shape_corr->SetXTitle("FADC sample no. + ( - t_{arrival}) / T_{ADC}"); + shape_corr->SetYTitle("signal [FADC counts]"); + shape_corr->SetMarkerStyle(20); + shape_corr->SetMarkerColor(4); + shape_corr->SetMarkerSize(0.3); + shape_corr->Draw(); + + + TCanvas * c_sig = new TCanvas("c_sig","c_sig",600,400); + c_sig->SetFillColor(0); + c_sig->SetBorderMode(0); + c_sig->SetGridx(); + c_sig->SetGridy(); + + hsig->SetStats(0); + TString title4 = "Noise Extracted Charge, MExtractFixedWindowPeakSearch, Pixel "; + title4+=ipix; + hsig->SetTitle(title4); + hsig->SetXTitle("reco charge [FADC counts]"); + hsig->SetYTitle("events / 2 FADC counts"); + hsig->SetLineColor(4); + hsig->SetLineWidth(2); + hsig->Draw(); + + TCanvas * c_time = new TCanvas("c_time","c_time",600,400); + c_time->SetFillColor(0); + c_time->SetBorderMode(0); + c_time->SetGridx(); + c_time->SetGridy(); + + htime->SetStats(0); + TString title5 = "Noise Arrival Time, MExtractTimeHighestIntegral, Pixel "; + title5+=ipix; + htime->SetTitle(title5); + htime->SetXTitle("arrival time [T_{ADC}]"); + htime->SetYTitle("events / 0.1 T_{ADC}"); + htime->SetLineColor(4); + htime->SetLineWidth(2); + htime->Draw(); + */ +} + + + + + + +