Index: /fact/tools/rootmacros/GetDrsOffset.C =================================================================== --- /fact/tools/rootmacros/GetDrsOffset.C (revision 12581) +++ /fact/tools/rootmacros/GetDrsOffset.C (revision 12581) @@ -0,0 +1,195 @@ +/* +GetDrsOffset.C - a ROOT macro calculating the DRS Offset calibration constants + + -- TESTING -- + + In case you have taken a file, which is good for calibrating the DRS Offset, + then you can calculate the offset using this macro. + It just calculates the mean value and the RMS for all DRS bins of all pixels + in the file, for the number of events you like. + + The output is two TH2Fs, which contain all requested data on the one hand and + on the other hand are good for having an overview. + + call it like this: + GetDrsOffset("path/20111111_001.fits.gz") + + GetDrsOffset("path/20111111_001.fits.gz", "foo/output.root") + saves the output to a root file + + GetDrsOffset("path/20111111_001.fits.gz", "", "bar/test.bin") + saves the output to a binary file. + + GetDrsOffset("path/20111111_001.fits.gz", "foo/out.root", "bar/test.bin", false, -1) + produces all output, but does not open a ROOT Canvas. + -1 : means process all events in the given file. + +*/ + + +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include +#include + +#define HAVE_ZLIB +#include "fits.h" + +#include "openFits.h" +#include "openFits.c" + + +int NEvents; +vector AllPixelDataVector; +vector StartCellVector; +unsigned int CurrentEventID; +size_t PXLxROI; +UInt_t RegionOfInterest; +UInt_t NumberOfPixels; + +TH2F * hOffset; +TH2F * hOffsetRMS; +double mean[1440][1024]={{0}}; +int n[1440][1024]={{0}}; +double rms[1440][1024]={{0}}; +TObjArray hList; + + +int GetDrsOffset( + char *datafilename = "data/20111016_013.fits.gz", + const char *RootFileName = "", + const char *OutFileName = "", + bool graphic = true, + int nevents = -1, + int verbosityLevel = 1 + ) +{ +gStyle->SetPalette(1,0); + + + fits * datafile; + NEvents = openDataFits( datafilename, &datafile, + AllPixelDataVector, StartCellVector, CurrentEventID, + RegionOfInterest, NumberOfPixels, PXLxROI, verbosityLevel); + if (NEvents == 0){ + cout << "return code of OpenDataFile:" << datafilename<< endl; + cout << "is zero -> aborting." << endl; + return 1; + } + + if (verbosityLevel > 0) + cout <<"number of events in file: "<< NEvents << "\t"; + if ( nevents == -1 || nevents > NEvents ) nevents = NEvents; // -1 means all! + if (verbosityLevel > 0) + cout <<"of, which "<< nevents << " will be processed"<< endl; + + int npixel = -1; + if (verbosityLevel > 0) + cout <<"Total # of Pixel: "<< NumberOfPixels << "\t"; + if ( npixel == -1 || npixel > (int)NumberOfPixels ) npixel = (int)NumberOfPixels; // -1 means all! + if (verbosityLevel > 0) + cout <<"of, which "<< npixel << " will be processed"<< endl; + + + if (verbosityLevel > 0) + cout <<" The ROI of this datafile is:" << RegionOfInterest << endl; + + double cur; + // -------------------- OFFSET -------------------------------------------- + for ( int ev = 0; ev < nevents; ev++) { + printf("%5d / %5d\b\b\b\b\b\b\b\b\b\b\b\b\b",ev,nevents); + datafile->GetRow( ev ); + for ( int pix = 0 ; pix < npixel; pix++ ){ + for (unsigned int sl=0; sl< RegionOfInterest; sl++){ + cur = AllPixelDataVector[pix * RegionOfInterest + sl]; + mean[pix][(sl+StartCellVector[pix])%1024] += cur; + ++n[pix][(sl+StartCellVector[pix])%1024]; + rms[pix][(sl+StartCellVector[pix])%1024] += cur * cur; + } + } + } + for (int p=0; p<1440; p++) + for (int s=0; s<1024; s++){ + if (n[p][s]==0){ + cout << "Error: no data for pixel:" << p <<" slice:" << s << endl; + return 1; + } + mean[p][s]/=n[p][s]; + rms[p][s]/=n[p][s]; + rms[p][s]-=mean[p][s]*mean[p][s]; + rms[p][s]=sqrt(rms[p][s]); + } + + + // ----------------- For plotting only ------------------- + if (graphic){ + + TCanvas *cOffset = NULL; + cOffset = new TCanvas("cOffset","Offset",10,10,800,400); + cOffset->Divide(1, 2); + + hOffset = new TH2F("hOffset","Offset",1440,-0.5,1439.5, 1024, -0.5, 1023.5); + hOffset->GetXaxis()->SetTitle( "pixel" ); + hOffset->GetYaxis()->SetTitle( "slice" ); + hList.Add( hOffset ); + + hOffsetRMS = new TH2F("hOffsetRMS","Offset RMS",1440,-0.5,1439.5, 1024, -0.5, 1023.5); + hOffsetRMS->GetXaxis()->SetTitle( "pixel" ); + hOffsetRMS->GetYaxis()->SetTitle( "slice" ); + hList.Add( hOffsetRMS ); + + for ( int pix = 0 ; pix < 1440; pix++ ){ + for (unsigned int sl=0; sl< 1024; sl++){ + hOffset->Fill(pix,sl,mean[pix][sl]); + hOffsetRMS->Fill(pix,sl,rms[pix][sl]); + } + } + cOffset->cd(1); + hOffset->Draw("COLZ"); + cOffset->cd(2); + hOffsetRMS->Draw("COLZ"); + cOffset->Modified(); + cOffset->Update(); + } + + if ( strlen(OutFileName) > 0 ) { + ofstream out(OutFileName, ios::out | ios::binary); + if(!out) { + cout << "Cannot open file." << OutFileName << endl; + return 1; + } + out.write((char *) &mean, sizeof mean); + out.write((char *) &rms, sizeof rms); + out.close(); + cout << "binary output written to \"" << OutFileName << "\"" << endl; + }else + cout << "no output generated" << endl; + +/* + double offset[1440][1024]; + ifstream in(OutFileName, ios::in | ios::binary); + in.read((char*) &offset, sizeof offset); + in.close(); + + + for ( int pix = 0 ; pix < 1440; pix++ ){ + for (unsigned int sl=0; sl< 1024; sl++){ + cout << offset[pix][sl] << endl; + } + } +*/ + + TFile tf( RootFileName, "RECREATE"); + hList.Write(); // write the major histograms into the top level directory + tf.Close(); // close the file + return( 0 ); +} + Index: /fact/tools/rootmacros/README.txt =================================================================== --- /fact/tools/rootmacros/README.txt (revision 12580) +++ /fact/tools/rootmacros/README.txt (revision 12581) @@ -151,3 +151,29 @@ +--------------------------------------------------------------------------------------- +GetDrsOffset.C - a ROOT macro calculating the DRS Offset calibration constants + + -- TESTING -- + In case you have taken a file, which is good for calibrating the DRS Offset, + then you can calculate the offset using this macro. + It just calculates the mean value and the RMS for all DRS bins of all pixels + in the file, for the number of events you like. + + The output is two TH2Fs, which contain all requested data on the one hand and + on the other hand are good for having an overview. + + call it like this: + GetDrsOffset("path/20111111_001.fits.gz") + + GetDrsOffset("path/20111111_001.fits.gz", "foo/output.root") + saves the output to a root file + + GetDrsOffset("path/20111111_001.fits.gz", "", "bar/test.bin") + saves the output to a binary file. + + GetDrsOffset("path/20111111_001.fits.gz", "foo/out.root", "bar/test.bin", false, -1) + produces all output, but does not open a ROOT Canvas. + -1 : means process all events in the given file. + +---------------------------------------------------------------------------------------