git-svn-id: file:///afs/psi.ch/project/sls_det_software/svn/slsDetectorCalibration@49 113b152e-814d-439b-b186-022a431db7b5

This commit is contained in:
jungmann_j 2014-04-08 14:17:49 +00:00
parent 9297b14500
commit 9b05e59ebc
12 changed files with 207 additions and 328 deletions

View File

@ -447,12 +447,10 @@ TF1* energyCalibration::fitFunction(TF1 *fun, TH1 *h1, Double_t *mypar, Double_t
fitfun= h1->GetFunction(fname); fitfun= h1->GetFunction(fname);
if (fitfun) {
fitfun->GetParameters(mypar); fitfun->GetParameters(mypar);
for (int ip=0; ip<6; ip++) { for (int ip=0; ip<6; ip++) {
emypar[ip]=fitfun->GetParError(ip); emypar[ip]=fitfun->GetParError(ip);
} }
}
return fitfun; return fitfun;
} }

View File

@ -1,19 +1,11 @@
void gMap(char *tit, float g=1) { {
//.L energyCalibration.cpp+ //.L energyCalibration.cpp+
//.L gainMap.C+ //.L gainMap.C+
TFile fin("/data/moench_xbox_20140113/MoTarget_45kV_0_8mA_120V_cds_g4.root");
char fname[1000];
sprintf(fname,"/data/moench_xbox_20140113/MoTarget_45kV_0_8mA_120V_%s.root",tit);
TFile fin(fname);
TH2F *h2=fin.Get("h2"); TH2F *h2=fin.Get("h2");
TH2F *gMap=(TH2F*)gainMap(h2,g); TH2F *gMap=gainMap(h2,4);
gMap->Draw("colz"); gMap->Draw("colz");
sprintf(fname,"/data/moench_xbox_20140113/gain_map_%s.root",tit);
TFile fout(fname,"RECREATE");
gMap->Write();
fout.Close();
} }

View File

@ -16,25 +16,22 @@ using namespace std;
#define FOPT "0" #define FOPT "0"
THStack *gainMap(TH2F *h2, float g) { TH2F *gainMap(TH2F *h2, float g) {
// int npx, npy; // int npx, npy;
int npx=160, npy=160; int npx=160, npy=160;
// det->getDetectorSize(npx, npy); // det->getDetectorSize(npx, npy);
THStack *hs=new THStack();
TH2F *gMap=new TH2F("gmap",h2->GetTitle(), npx, -0.5, npx-0.5, npy, -0.5, npy-0.5); TH2F *gMap=new TH2F("gmap",h2->GetTitle(), npx, -0.5, npx-0.5, npy, -0.5, npy-0.5);
TH2F *nMap=new TH2F("nmap",h2->GetTitle(), npx, -0.5, npx-0.5, npy, -0.5, npy-0.5);
hs->Add(gMap); Double_t ens[3]={0,8,17.5}, eens[3]={0.,0.1,0.1};
hs->Add(nMap); Double_t peaks[3], epeaks[3];
Double_t ens[1]={20.}, eens[1]={20.};
Double_t peaks[1], epeaks[1];
int ibin; int ibin;
TH1D *px; TH1D *px;
energyCalibration *enCal=new energyCalibration(); energyCalibration *enCal=new energyCalibration();
enCal->setPlotFlag(0); enCal->setPlotFlag(0);
// enCal->setChargeSharing(0); // enCal->setChargeSharing(0);
@ -53,44 +50,90 @@ THStack *gainMap(TH2F *h2, float g) {
TGraph *glin; TGraph *glin;
Double_t peakdum, hpeakdum; Double_t peakdum, hpeakdum;
int ix=20; for (int ix=1; ix<npx-1; ix++) {
int iy=40; for (int iy=1; iy<npy-1; iy++) {
for ( ix=1; ix<npx-1; ix++) {
for ( iy=1; iy<npy-1; iy++) {
// cout << ix << " " << iy << " " << ibin << endl; // cout << ix << " " << iy << " " << ibin << endl;
ibin=ix*npy+iy; ibin=ix*npy+iy;
px=h2->ProjectionX("px",ibin+1,ibin+1); px=h2->ProjectionX("px",ibin+1,ibin+1);
enCal->setFitRange(50,3000); prms=10*g;
//enCal->setChargeSharing(0); iit=0;
if (px) { np=s->Search(px,prms,"",0.2);
enCal->fixParameter(0,0); while (np !=2) {
enCal->fixParameter(1,0); if (np>2)
enCal->fixParameter(5,0); prms+=0.5*prms;
else
prms-=0.5*prms;
iit++;
if (iit>=10)
break;
np=s->Search(px,prms,"",0.2);
}
if (np!=2)
cout << "peak search could not converge " << ibin << endl;
if (np==2) {
pm=NULL;
functions=px->GetListOfFunctions();
if (functions)
pm = (TPolyMarker*)functions->FindObject("TPolyMarker");
if (pm) {
peakX=pm->GetX();
peakY=pm->GetY();
if (peakX[0]>peakX[1]) {
peakdum=peakX[0];
hpeakdum=peakY[0];
peakX[0]= peakX[1];
peakY[0]= peakY[1];
peakX[1]= peakdum;
peakY[1]= hpeakdum;
}
cout << "("<< ix << "," << iy << ") " << endl;
for (int ip=0; ip<np; ip++) {
// cout << ip << " " << peakX[ip] << " " << peakY[ip] << endl;
enCal->setFitRange(peakX[ip]-10*g,peakX[ip]+10*g);
mypar[0]=0; mypar[0]=0;
mypar[1]=0; mypar[1]=0;
mypar[2]=px->GetBinCenter(px->GetMaximumBin()); mypar[2]=peakX[ip];
mypar[3]=10; mypar[3]=g*10;
mypar[4]=px->GetMaximum(); mypar[4]=peakY[ip];
mypar[5]=0; mypar[5]=0;
enCal->setStartParameters(mypar); enCal->setStartParameters(mypar);
enCal->fitSpectrum(px,mypar,emypar); enCal->fitSpectrum(px,mypar,emypar);
cout << ix << " " << iy << " " << mypar[2] << endl;
peaks[ip+1]=mypar[2];
epeaks[ip+1]=emypar[2];
} }
if (mypar[2]>0) { peaks[0]=0;
gMap->SetBinContent(ix+1,iy+1,mypar[2]/ens[0]); epeaks[0]=1;
gMap->SetBinError(ix+1,iy+1,emypar[2]/ens[0]);
nMap->SetBinContent(ix+1,iy+1,mypar[3]); // for (int i=0; i<3; i++) cout << i << " " << ens[i] << " " << eens[i]<< " "<< peaks[i]<< " " << epeaks[i] << endl;
nMap->SetBinError(ix+1,iy+1,emypar[3]);
glin= enCal->linearCalibration(3,ens,eens,peaks,epeaks,gain,off,egain,eoff);
// cout << "Gain " << gain << " off " << off << endl;
if (off>-10 && off<10) {
gMap->SetBinContent(ix+1,iy+1,gain);
gMap->SetBinError(ix+1,iy+1,egain);
} }
if (glin)
delete glin;
} }
} }
return hs;
}
}
return gMap;
} }
@ -107,48 +150,56 @@ TH2F *noiseMap(TH2F *h2) {
for (int ix=0; ix<npx; ix++) { for (int ix=0; ix<npx; ix++) {
for (int iy=0; iy<npy; iy++) { for (int iy=0; iy<npy; iy++) {
cout << ix << " " << iy << " " << ibin << endl; cout << ix << " " << iy << " " << ibin << endl;
ibin=h2->GetYaxis()->FindBin(ix*npy+iy); ibin=ix*npy+iy;
px=h2->ProjectionX("px",ibin,ibin); px=h2->ProjectionX("px",ibin+1,ibin+1);
px->Fit("gaus", FOPT,"",-100,100); px->Fit("gaus", FOPT);
if (px->GetFunction("gaus")) { if (px->GetFunction("gaus")) {
if (px->GetFunction("gaus")->GetParameter(1)>-5 && px->GetFunction("gaus")->GetParameter(1)<5)
nMap->SetBinContent(ix+1,iy+1,px->GetFunction("gaus")->GetParameter(2)); nMap->SetBinContent(ix+1,iy+1,px->GetFunction("gaus")->GetParameter(2));
} }
// delete px; // delete px;
} }
} }
return nMap; return nMap;
} }
THStack *noiseHistos(TH2F *nmap, TH2F *gmap=NULL) {
THStack *noiseHistos(char *tit) {
char fname[10000];
char tit[1000]; sprintf(fname,"/data/moench_xbox_20140116/noise_map_%s.root",tit);
TFile *fn=new TFile(fname);
TH2F *nmap=(TH2F*)fn->Get("nmap");
if (nmap==NULL) { if (nmap==NULL) {
cout << "No noise map" << endl; cout << "No noise map in file " << fname << endl;
return NULL; return NULL;
} }
if (gmap) { sprintf(fname,"/data/moench_xbox_20140113/gain_map_%s.root",tit);
nmap->Divide(gmap); TFile *fg=new TFile(fname);
nmap->Scale(1000./3.6); TH2F *gmap=(TH2F*)fg->Get("gmap");
if (gmap==NULL) {
cout << "No gain map in file " << fname << endl;
return NULL;
} }
nmap->Divide(gmap);
nmap->Scale(1000./3.6);
strcpy(tit,nmap->GetTitle());
THStack *hs=new THStack(tit,tit); THStack *hs=new THStack(tit,tit);
hs->SetTitle(tit); hs->SetTitle(tit);
TH1F *h; TH1F *h;
char hname[100]; char hname[100];
cout << tit << endl; cout << tit << endl;
for (int is=0; is<4; is++) { for (int is=0; is<4; is++) {
sprintf(hname,"h%ds",is+1); sprintf(hname,"h%ds",is+1);
h=new TH1F(hname,tit,500,0,500); h=new TH1F(hname,tit,500,0,500);
hs->Add(h); hs->Add(h);
// cout << hs->GetHists()->GetEntries() << endl; // cout << hs->GetHists()->GetEntries() << endl;
@ -164,7 +215,7 @@ THStack *noiseHistos(TH2F *nmap, TH2F *gmap=NULL) {
h->SetLineColor(is+1); h->SetLineColor(is+1);
if (h->GetFunction("gaus")) { if (h->GetFunction("gaus")) {
h->GetFunction("gaus")->SetLineColor(is+1); h->GetFunction("gaus")->SetLineColor(is+1);
cout << " or " << h->GetFunction("gaus")->GetParameter(1) << "+-" << h->GetFunction("gaus")->GetParameter(2); cout << " or " << h->GetFunction("gaus")->GetParameter(1) << "+-" << h->GetFunction("gaus")->GetParError(1);
} }
cout << endl; cout << endl;
} }

View File

@ -20,7 +20,6 @@
//#include "MovingStat.h" //#include "MovingStat.h"
using namespace std; using namespace std;
#define NC 48 #define NC 48
@ -51,15 +50,14 @@ Loops over data file to find single photons, fills the tree (and writes it to fi
\param ymax maximum y coordinate \param ymax maximum y coordinate
\param cmsub enable commonmode subtraction \param cmsub enable commonmode subtraction
\param hitfinder if 0: performs pedestal subtraction, not hit finding; if 1: performs both pedestal subtraction and hit finding \param hitfinder if 0: performs pedestal subtraction, not hit finding; if 1: performs both pedestal subtraction and hit finding
\param fcal -> calibration file, which contains 7 columns per pixel; 1) pixel #, 2) cal. offset G0, 3) cal. slope G0, 4) cal. offset G1, 5) cal. slope G1, 6) cal. offset G2, 7) cal. slope G2
\returns pointer to histo stack with cluster spectra \returns pointer to histo stack with cluster spectra
*/ */
THStack *jungfrauReadData(char *fformat, char *tit, int runmin, int runmax, int nbins=1500, int hmin=-500, int hmax=1000, int sign=1, double hc=0, int xmin=1, int xmax=NC-1, int ymin=1, int ymax=NR-1, int cmsub=0, int hitfinder=1){ THStack *jungfrauReadData(char *fformat, char *tit, int runmin, int runmax, int nbins=1500, int hmin=-500, int hmax=1000, int sign=1, double hc=0, int xmin=1, int xmax=NC-1, int ymin=1, int ymax=NR-1, int cmsub=0, int hitfinder=1){
/*
// read in calibration file // read in calibration file
ifstream calfile("/mnt/slitnas/datadir_jungfrau02/analysis_tests/CalibrationParametersTest_fake.txt"); ifstream calfile("/home/l_msdetect/TriesteBeam2014/dummy data for scripts/CalibrationParametersTest_fake.txt");
if (calfile.is_open()==0){cout << "Unable to open calibration file!" << endl;} if (calfile.is_open()==0){cout << "Unable to open calibration file!" << endl;}
int pix; int pix;
double of0,sl0,of1,sl1,of2,sl2; double of0,sl0,of1,sl1,of2,sl2;
@ -73,10 +71,10 @@ THStack *jungfrauReadData(char *fformat, char *tit, int runmin, int runmax, int
sl_2[pix]=sl2; //if(pix==200) cout << "sl_2[200] " << sl_2[200] << endl; sl_2[pix]=sl2; //if(pix==200) cout << "sl_2[200] " << sl_2[200] << endl;
} }
calfile.close(); calfile.close();
*/
double adc_value, num_photon; double adc_value, num_photon;
jungfrau02Data *decoder=new jungfrau02Data(1,0,0);//(3,0,0); // (adc,offset,crosstalk) //(1,0,0) //(3,0,0) for readout of GB jungfrau02Data *decoder=new jungfrau02Data(1,0,0);//(1,0,0); // (adc,offset,crosstalk) //(1,0,0) //(3,0,0) for readout of GB
jungfrau02CommonMode *cmSub=NULL; jungfrau02CommonMode *cmSub=NULL;
if (cmsub) if (cmsub)
cmSub=new jungfrau02CommonMode(); cmSub=new jungfrau02CommonMode();
@ -107,9 +105,6 @@ THStack *jungfrauReadData(char *fformat, char *tit, int runmin, int runmax, int
TH2F *h1=new TH2F("h1",tit,nbins,hmin-0.5,hmax-0.5,NC*NR,-0.5,NC*NR-0.5); TH2F *h1=new TH2F("h1",tit,nbins,hmin-0.5,hmax-0.5,NC*NR,-0.5,NC*NR-0.5);
hs->Add(h1); hs->Add(h1);
if (hitfinder) { if (hitfinder) {
h2=new TH2F("h2",tit,nbins,hmin-0.5,hmax-0.5,NC*NR,-0.5,NC*NR-0.5); h2=new TH2F("h2",tit,nbins,hmin-0.5,hmax-0.5,NC*NR,-0.5,NC*NR-0.5);
h3=new TH2F("h3",tit,nbins,hmin-0.5,hmax-0.5,NC*NR,-0.5,NC*NR-0.5); h3=new TH2F("h3",tit,nbins,hmin-0.5,hmax-0.5,NC*NR,-0.5,NC*NR-0.5);
@ -193,7 +188,7 @@ THStack *jungfrauReadData(char *fformat, char *tit, int runmin, int runmax, int
thisEvent=filter->getEventType(buff, ix, iy, cmsub); thisEvent=filter->getEventType(buff, ix, iy, cmsub);
int gainBits=decoder->getGainBits(buff,ix,iy); // get gain bits int gainBits=decoder->getGainBits(buff,ix,iy); //XXX
#ifdef MY_DEBUG #ifdef MY_DEBUG
if (hitfinder) { if (hitfinder) {
@ -211,7 +206,6 @@ THStack *jungfrauReadData(char *fformat, char *tit, int runmin, int runmax, int
h1->Fill(filter->getClusterTotal(1), iy+NR*ix); h1->Fill(filter->getClusterTotal(1), iy+NR*ix);
if (hitfinder) { if (hitfinder) {
if (thisEvent==PHOTON_MAX ) { if (thisEvent==PHOTON_MAX ) {

View File

@ -66,6 +66,9 @@ class moench02ModuleData : public slsReceiverData<uint16_t> {
setDataMap(dMap); setDataMap(dMap);
setDataMask(dMask); setDataMask(dMask);
}; };
@ -129,45 +132,6 @@ class moench02ModuleData : public slsReceiverData<uint16_t> {
}; };
class moench02OversampledData : public moench02ModuleData {
public:
moench02OversampledData(int nos, int off=0): moench02ModuleData(), nSamples(nos), offset(off){};
/**
returns the pixel value as double correcting for the output buffer crosstalk
\param data pointer to the memory
\param ix coordinate in the x direction
\param iy coordinate in the y direction
\returns channel value as double
*/
double getValue(char *data, int ix, int iy=0) {
double v=0, is=0;
int ix1, iy1, isc=ix/40, ip=ix%40, ip1=iy*nSamples*40+ix;
if (iy*nSamples<160) {
for (int i=offset; i<nSamples; i++) {
ip1=iy*nSamples*40+ip*nSamples+i-1;
iy1=ip1/40;
ix1=ip1%40+isc*40;
if (ix1>=0 && iy1>=0 && ix1<160 && iy1<160) {
v+=moench02ModuleData::getValue(data,ix1,iy1);
is++;
}
}
if (is>0)
v/=is;
}
return v;
};
private:
int nSamples;
int offset;
};
#endif #endif

View File

@ -79,7 +79,6 @@ THStack *moenchReadData(char *fformat, char *tit, int runmin, int runmax, int nb
TH2F *h3; TH2F *h3;
TH2F *hetaX; TH2F *hetaX;
TH2F *hetaY; TH2F *hetaY;
TH2D *clustHist;
THStack *hs=new THStack("hs",fformat); THStack *hs=new THStack("hs",fformat);
@ -103,24 +102,19 @@ THStack *moenchReadData(char *fformat, char *tit, int runmin, int runmax, int nb
ifstream filebin; ifstream filebin;
int iev;
int ix=20, iy=20, ir, ic; int ix=20, iy=20, ir, ic;
Int_t iFrame; Int_t iFrame;
TTree *tall; TTree *tall;
cout << "init tree " << tit << endl;
if (hitfinder) if (hitfinder)
tall=filter->initEventTree(tit, &iFrame); tall=filter->initEventTree(tit, &iFrame);
cout << "done" << endl;
#ifdef MY_DEBUG #ifdef MY_DEBUG
quadrant quad;
tall->Branch("q",&quad,"q/I");
TCanvas *myC; TCanvas *myC;
TH2F *he; TH2F *he;
@ -128,9 +122,6 @@ THStack *moenchReadData(char *fformat, char *tit, int runmin, int runmax, int nb
TCanvas *cH2; TCanvas *cH2;
TCanvas *cH3; TCanvas *cH3;
int quadrants[5];
for(int i = 0; i < 5; i++){ quadrants[i] = 0; }
if (hitfinder) { if (hitfinder) {
myC=new TCanvas(); myC=new TCanvas();
he=new TH2F("he","Event Mask",xmax-xmin, xmin, xmax, ymax-ymin, ymin, ymax); he=new TH2F("he","Event Mask",xmax-xmin, xmin, xmax, ymax-ymin, ymin, ymax);
@ -145,8 +136,6 @@ THStack *moenchReadData(char *fformat, char *tit, int runmin, int runmax, int nb
cH3=new TCanvas(); cH3=new TCanvas();
cH3->SetLogz(); cH3->SetLogz();
h3->Draw("colz"); h3->Draw("colz");
clustHist= new TH2D("clustHist","clustHist",3,-1.5,1.5,3,-1.5,1.5);
} }
#endif #endif
filter->newDataSet(); filter->newDataSet();
@ -154,7 +143,7 @@ THStack *moenchReadData(char *fformat, char *tit, int runmin, int runmax, int nb
for (int irun=runmin; irun<runmax; irun++) { for (int irun=runmin; irun<runmax; irun++) {
sprintf(fname,fformat,irun); sprintf(fname,fformat,irun);
cout << "file name " << fname << " ( " << (((double)(irun-runmin))*100./((double)(runmax-runmin))) << "% )" << endl; cout << "file name " << fname << endl;
filebin.open((const char *)(fname), ios::in | ios::binary); filebin.open((const char *)(fname), ios::in | ios::binary);
nph=0; nph=0;
while ((buff=decoder->readNextFrame(filebin))) { while ((buff=decoder->readNextFrame(filebin))) {
@ -198,17 +187,9 @@ THStack *moenchReadData(char *fformat, char *tit, int runmin, int runmax, int nb
h3->Fill(filter->getClusterTotal(3), iy+NR*ix); h3->Fill(filter->getClusterTotal(3), iy+NR*ix);
iFrame=decoder->getFrameNumber(buff); iFrame=decoder->getFrameNumber(buff);
#ifdef MY_DEBUG
for(int cx=-1; cx <2; cx++)
for(int cy=-1; cy < 2; cy++)
clustHist->Fill(cx,cy,filter->getClusterElement(cx,cy));
quad = filter->getQuadrant();
quadrants[quad]++;
#endif
tall->Fill(); tall->Fill();
} }
@ -234,10 +215,10 @@ THStack *moenchReadData(char *fformat, char *tit, int runmin, int runmax, int nb
#endif #endif
nf++; nf++;
//cout << "=" ; cout << "=" ;
delete [] buff; delete [] buff;
} }
//cout << nph << endl; cout << nph << endl;
if (filebin.is_open()) if (filebin.is_open())
filebin.close(); filebin.close();
else else
@ -250,14 +231,6 @@ THStack *moenchReadData(char *fformat, char *tit, int runmin, int runmax, int nb
delete decoder; delete decoder;
cout << "Read " << nf << " frames" << endl; cout << "Read " << nf << " frames" << endl;
#ifdef MY_DEBUG
cout << "quadrants: " ;
for(int i = 0; i<5; i++) cout << i << ": " << quadrants[i] << " || " ;
cout << endl;
cout << "Read Events " << nph << endl;
#endif
return hs; return hs;
} }

View File

@ -9,30 +9,23 @@ typedef struct task_s{
int runmin; int runmin;
int runmax; int runmax;
int treeIndex; int treeIndex;
int tNumber;
double xTalkFactor;
} Task; } Task;
double hc = 0.1; // read - out crosstalk
void *moenchMakeTreeTask(void *p){ void *moenchMakeTreeTask(void *p){
TThread::Lock(); TThread::Lock();
char fname[1000]; char fname[1000];
Task *t = (Task *)p; Task *t = (Task *)p;
sprintf(fname,"%s%s_%i.root",t->tdir,t->tname,t->treeIndex); sprintf(fname,"%s%s_%i.root",t->tdir,t->tname,t->treeIndex);
TFile *f = new TFile(fname,"RECREATE"); TFile *f = new TFile(fname,"RECREATE");
double xTalkC = t->xTalkFactor; //((double)t->tNumber) * 0.01; Calibrated with data in /data/moench_trieste_20140313_trees/xtalkScan cout << "Call moenchReadData(" << t->fformat << "," << t->tname << "," << t->runmin<< "," << t->runmax <<")" << endl;
cout << "Call moenchReadData(" << t->fformat << "," << t->tname << "," << t->runmin<< "," << t->runmax <<") with xTalkC: " << xTalkC << endl;
TThread::UnLock(); TThread::UnLock();
moenchReadData(t->fformat,t->tname,t->runmin,t->runmax, 1500, -500, 1000, 1, xTalkC); moenchReadData(t->fformat,t->tname,t->runmin,t->runmax);
if(f && f->IsOpen()){
f->Close(); f->Close();
}
return 0; return 0;
} }
void moenchReadDataMT(char *fformat, char *tit, char *tdir, int runmin, int runoffset, int nThreads, int treeIndexStart=0, double xTalkFactor=0.044){ void moenchReadDataMT(char *fformat, char *tit, char *tdir, int runmin, int runoffset, int nThreads, int treeIndexStart=0){
char threadName[1000]; char threadName[1000];
TThread *threads[nThreads]; TThread *threads[nThreads];
for(int i = 0; i < nThreads; i++){ for(int i = 0; i < nThreads; i++){
@ -44,8 +37,6 @@ void moenchReadDataMT(char *fformat, char *tit, char *tdir, int runmin, int runo
t->runmin = runmin + i*runoffset; t->runmin = runmin + i*runoffset;
t->runmax = runmin + (i+1)*runoffset - 1; t->runmax = runmin + (i+1)*runoffset - 1;
t->treeIndex = treeIndexStart + i; t->treeIndex = treeIndexStart + i;
t->tNumber = i;
t->xTalkFactor = xTalkFactor;
cout << "start thread " << i << " start: " << t->runmin << " end " << t->runmax << endl; cout << "start thread " << i << " start: " << t->runmin << " end " << t->runmax << endl;
threads[i] = new TThread(threadName, moenchMakeTreeTask, t); threads[i] = new TThread(threadName, moenchMakeTreeTask, t);
threads[i]->Run(); threads[i]->Run();
@ -55,32 +46,7 @@ void moenchReadDataMT(char *fformat, char *tit, char *tdir, int runmin, int runo
threads[i]->Join(); threads[i]->Join();
} }
cout << " ( DONE ) " << endl;
} }
//to compile: g++ -DMYROOT -DMYROOT1 -g `root-config --cflags --glibs` -o moenchReadDataMT moenchReadDataMT.C
int main(int argc, char **argv){
if(argc < 8){
cout << "Usage: " << argv[0] << " fformat tit tdir runmin runoffset nThreads treeIndexStart [xTalkFactor]" << endl;
exit(-1);
}
char *fformat = argv[1];
char *tit = argv[2];
char *tdir = argv[3];
int runmin = atoi(argv[4]);
int runoffset = atoi(argv[5]);
int nThreads = atoi(argv[6]);
int treeIndexStart = atoi(argv[7]);
double xTalkFactor = 0.044;
if(argc == 9)
xTalkFactor = atof(argv[8]);
moenchReadDataMT(fformat, tit, tdir,runmin,runoffset,nThreads,treeIndexStart, xTalkFactor);
}

View File

@ -1,5 +1,4 @@
#include "moenchReadData.C" #include "moenchReadData.C"
#include "moenchReadOversampledData.C"
@ -14,12 +13,12 @@ void raedNoiseData(char *tit, int sign=1){
TFile *fout; TFile *fout;
THStack *hs2N; THStack *hs2N;
sprintf(fname,"/data/moench_trieste_calibration_trees/flat_20keV_%s.root",tit); sprintf(fname,"/data/moench_xbox_20140113/MoTarget_45kV_0_8mA_120V_%s_0.root",tit);
fout=new TFile(fname,"RECREATE"); fout=new TFile(fname,"RECREATE");
sprintf(fname,"/data/moench_trieste_calibration/flat_20keV_%s_f00000%%04d000_0.raw",tit); sprintf(fname,"/data/moench_xbox_20140113/MoTarget_45kV_0_8mA_120V_%s_f00000%%04d000_0.raw",tit);
hs2N=moenchReadData(fname,tit,0,2000,1500,-500,2500,sign,0.,1,159,1,159, 0,1); hs2N=moenchReadData(fname,0,3000,1500,-500,2500,1,0.,1,159,1,159, 0,1);
hs2N->SetName(tit); hs2N->SetName(tit);
hs2N->SetTitle(tit); hs2N->SetTitle(tit);
(TH2F*)(hs2N->GetHists()->At(0))->Write(); (TH2F*)(hs2N->GetHists()->At(0))->Write();
@ -37,59 +36,6 @@ void raedNoiseData(char *tit, int sign=1){
} }
void raedOsNoiseData(){
THStack *hs[10];
char fname[1000];
char f[1000];
TFile *fout;
strcpy(fname,"/data/moench_xbox_201401_trees/noise_cds_g1_os10_2ndsample.root");
fout=new TFile(fname,"RECREATE");
for (int ir=0; ir<10; ir++) {
sprintf(fname,"/data/moench_xbox_201401/moench_xbox_20140116/noise_os10_16rows_f00000%%04d000_%d.raw",ir);
hs[ir]=moenchReadOversampledData(fname,"cds_g1_os10",0,1000,1000,-500,500,1,10,0,160,0,16,0);
}
TH2F *h;
int ii=0, ii1=0;
h=(TH2F*)(hs[0]->GetHists()->At(0));
TH2F *hn=(TH2F*)h->Clone();
for (int ir=1; ir<10; ir++) {
h=(TH2F*)(hs[ir]->GetHists()->At(0));
for (int iy=0; iy<16; iy++)
for (int ix=0; ix<160; ix++)
for (int ib=0; ib<hn->GetNbinsX(); ib++) {
ii=h->GetYaxis()->FindBin(iy+ix*160);
ii1= hn->GetYaxis()->FindBin(iy+ir*16+ix*160);
hn->SetBinContent(ib+1,ii1,h->GetBinContent(ib+1,ii));
}
}
hn->Write();
// (TH2F*)(hs2N->GetHists()->At(1))->Write();
// (TH2F*)(hs2N->GetHists()->At(2))->Write();
// (TH2F*)(hs2N->GetHists()->At(3))->Write();
// (TH2F*)(hs2N->GetHists()->At(4))->Write();
fout->Close();
}
void raedNoiseDataN(char *tit, int sign=1){ void raedNoiseDataN(char *tit, int sign=1){
@ -100,10 +46,10 @@ void raedNoiseDataN(char *tit, int sign=1){
TFile *fout; TFile *fout;
THStack *hs2N; THStack *hs2N;
sprintf(fname,"/data/moench_noise_20140327_trees/noise_%s.root",tit); sprintf(fname,"/data/moench_xbox_20140116/noise_%s.root",tit);
fout=new TFile(fname,"RECREATE"); fout=new TFile(fname,"RECREATE");
sprintf(fname,"/data/moench_noise_20140327/noise_%s_f00000%%04d000_0.raw",tit); sprintf(fname,"/data/moench_xbox_20140116/noise_%s_f00000%%04d000_0.raw",tit);
hs2N=moenchReadData(fname,tit,0,3000,1500,-500,2500,sign,0.,1,159,1,159, 0,0); hs2N=moenchReadData(fname,tit,0,3000,1500,-500,2500,sign,0.,1,159,1,159, 0,0);
hs2N->SetName(tit); hs2N->SetName(tit);
@ -128,7 +74,7 @@ void g4() {
raedNoiseData("cds_g4_low_gain"); raedNoiseData("cds_g4_low_gain");
raedNoiseData("cds_g4_sto1_only"); raedNoiseData("cds_g4_sto1_only");
raedNoiseData("cds_g4_no_sto"); raedNoiseData("cds_g4_no sto");
@ -138,7 +84,7 @@ void no_cds() {
raedNoiseData("cds_disable_low_gain",-1); raedNoiseData("cds_disable_low_gain",-1);
raedNoiseData("cds_disable_sto1_only",-1); raedNoiseData("cds_disable_sto1_only",-1);
raedNoiseData("cds_disable_no_sto",-1); raedNoiseData("cds_disable_no sto",-1);

View File

@ -5,14 +5,14 @@
TFile *fout; // output file TFile *fout; // output file
THStack *hs2N; // collection of objects THStack *hs2N; // collection of objects
fout=new TFile("/mnt/slitnas/datadir_jungfrau02/jungfrau_Trieste_2014/TriesteBeam2014/Analysis/D4_M1_HG_2000clk_CD986_CMsub.root","RECREATE"); // //fout=new TFile("/mnt/slitnas/datadir_jungfrau02/analysis_tests/tests_random.root","RECREATE"); //
fout=new TFile("/mnt/slitnas/datadir_jungfrau02/20140404_lowEXrays/test.root","RECREATE");
//hs2N=jungfrauReadData("/mnt/slitnas/datadir_jungfrau02/20140131_BeamSLS/datadir_310114/21keV_10us_HIGHG0_CDS550_%d.bin","tit",0,75,1500,-500,2500,1,0.,1,47,1,47,1,1); //1,1); hs2N=jungfrauReadData("/mnt/slitnas/datadir_jungfrau02/20140404_lowEXrays/Ti_1000clk_13kV40mA_%d.bin","tit",0,86,3000,-499.5,2500.5,1,0,1,47,1,47,1,1); //1,1);
//hs2N=jungfrauReadData("/mnt/slitnas/datadir_jungfrau02/20140228_GainBits/GSdata_laser_%d.bin","tit",0,5,1500,-500,2500,1,0.,1,47,1,47,0,1);//,"/mnt/slitnas/datadir_jungfrau02/analysis_tests/CalibrationParametersTest_fake.txt"); //1,1); // data set test GB -> set (3,0,0) in junfrauReadData.C //hs2N=jungfrauReadData("/mnt/slitnas/datadir_jungfrau02/20140228_GainBits/GSdata_laser_%d.bin","tit",0,5,1500,-500,2500,1,0.,1,47,1,47,0,1);//,"/mnt/slitnas/datadir_jungfrau02/analysis_tests/CalibrationParametersTest_fake.txt"); //1,1); // data set test GB -> set (3,0,0) in junfrauReadData.C
hs2N=jungfrauReadData("/mnt/slitnas/datadir_jungfrau02/jungfrau_Trieste_2014/Trieste_20140314/HG_20keV_CDS986_2000clk_%d.bin","tit",0,201,3000,-499.5,2500.5,1,0.,1,47,1,47,1,1); //hs2N=jungfrauReadData("/mnt/slitnas/datadir_jungfrau02/digital_bit_adc_test/vb_1825_%d.bin","tit",0,7,1500,-500,2500,1,0.,1,47,1,47,0,1);//,"/mnt/slitnas/datadir_jungfrau02/analysis_tests/CalibrationParametersTest_fake.txt"); //1,1); // data set test GB -> set (3,0,0) in junfrauReadData.C
//,"/mnt/slitnas/datadir_jungfrau02/analysis_tests/CalibrationParametersTest_fake.txt"); //1,1); // data set test GB -> set (3,0,0) in junfrauReadData.C
//hs2N->SetName("cds_g4"); //hs2N->SetName("cds_g4");
//hs2N->SetTitle("cds_g4"); //hs2N->SetTitle("cds_g4");
@ -25,7 +25,7 @@
//(TH2F*)(hs2N->GetHists()->At(5))->Write(); //(TH2F*)(hs2N->GetHists()->At(5))->Write();
//(TH2F*)(hs2N->GetHists()->At(6))->Write(); //(TH2F*)(hs2N->GetHists()->At(6))->Write();
//fout->Close(); // uncomment fout->Close(); // uncomment
} }

View File

@ -92,8 +92,6 @@ class singlePhotonDetector {
setClusterSize(csize); setClusterSize(csize);
}; };
/** /**
destructor. Deletes the cluster structure and the pdestalSubtraction array destructor. Deletes the cluster structure and the pdestalSubtraction array
*/ */
@ -335,27 +333,22 @@ class singlePhotonDetector {
\returns returns pointer to the TTree \returns returns pointer to the TTree
*/ */
TTree *initEventTree(char *tname, int *iFrame=NULL) { TTree *initEventTree(char *tname, int *iFrame=NULL) {
cout << tname << endl;
TTree* tall=new TTree(tname,tname); TTree* tall=new TTree(tname,tname);
cout << "tree instantiated" << endl;
if (iFrame) if (iFrame)
tall->Branch("iFrame",iFrame,"iframe/I"); tall->Branch("iFrame",iFrame,"iframe/I");
else else
tall->Branch("iFrame",&(cluster->iframe),"iframe/I"); tall->Branch("iFrame",&(cluster->iframe),"iframe/I");
cout << "iframe" << endl;
char tit[100];
tall->Branch("x",&(cluster->x),"x/I"); tall->Branch("x",&(cluster->x),"x/I");
tall->Branch("y",&(cluster->y),"y/I"); tall->Branch("y",&(cluster->y),"y/I");
char tit[100];
sprintf(tit,"data[%d]/D",clusterSize*clusterSizeY); sprintf(tit,"data[%d]/D",clusterSize*clusterSizeY);
tall->Branch("data",cluster->data,tit); tall->Branch("data",cluster->data,tit);
tall->Branch("pedestal",&(cluster->ped),"pedestal/D"); // tall->Branch("pedestal",&(cluster->ped),"pedestal/D");
tall->Branch("rms",&(cluster->rms),"rms/D"); // tall->Branch("rms",&(cluster->rms),"rms/D");
cout << "Cluster" << endl;
return tall; return tall;
}; };
#else
/** write cluster to filer*/
void writeCluster(FILE* myFile){cluster->write(myFile);};
#endif #endif
@ -380,6 +373,7 @@ class singlePhotonDetector {
double tot; /**< sum of the 3x3 cluster */ double tot; /**< sum of the 3x3 cluster */
double quadTot; /**< sum of the maximum 2x2cluster */ double quadTot; /**< sum of the maximum 2x2cluster */
}; };

View File

@ -193,7 +193,6 @@ class slsDetectorData {
*/ */
virtual double getValue(char *data, int ix, int iy=0) {return (double)getChannel(data, ix, iy);}; virtual double getValue(char *data, int ix, int iy=0) {return (double)getChannel(data, ix, iy);};
/** /**
Returns the frame number for the given dataset. Purely virtual func. Returns the frame number for the given dataset. Purely virtual func.
@ -205,15 +204,16 @@ class slsDetectorData {
virtual int getFrameNumber(char *buff)=0; virtual int getFrameNumber(char *buff)=0;
/* /**
Returns the packet number for the given dataset. purely virtual func Returns the packet number for the given dataset. purely virtual func
\param buff pointer to the dataset \param buff pointer to the dataset
\returns packet number number \returns packet number number
*/
virtual int getPacketNumber(char *buff)=0; virtual int getPacketNumber(char *buff)=0;
*/
/** /**

View File

@ -26,6 +26,7 @@ class slsReceiverData : public slsDetectorData<dataType> {
slsReceiverData(int npx, int npy, int np, int psize, int **dMap=NULL, dataType **dMask=NULL, int **dROI=NULL): slsDetectorData<dataType>(npx, npy, np*psize, dMap, dMask, dROI), nPackets(np), packetSize(psize) {}; slsReceiverData(int npx, int npy, int np, int psize, int **dMap=NULL, dataType **dMask=NULL, int **dROI=NULL): slsDetectorData<dataType>(npx, npy, np*psize, dMap, dMask, dROI), nPackets(np), packetSize(psize) {};
/** /**
Returns the frame number for the given dataset. Virtual func: works for slsDetectorReceiver data (also for each packet), but can be overloaded. Returns the frame number for the given dataset. Virtual func: works for slsDetectorReceiver data (also for each packet), but can be overloaded.
@ -44,7 +45,7 @@ class slsReceiverData : public slsDetectorData<dataType> {
*/ */
virtual int getPacketNumber(char *buff){return (*(int*)buff)&0xff;}; virtual int getPacketNumber(char *buff) {return (*(int*)buff)&0xff;};
@ -73,11 +74,11 @@ class slsReceiverData : public slsDetectorData<dataType> {
retval=NULL; retval=NULL;
np=0; np=0;
} else if (pnum==1) { } else if (pnum==1) {
fnum=fn;
retval=p; retval=p;
if (np>0) if (np>0)
/*cout << "*Incomplete frame number " << fnum << endl;*/ /*cout << "*Incomplete frame number " << fnum << endl;*/
np=0; np=0;
fnum=fn;
} else if (fn!=fnum) { } else if (fn!=fnum) {
if (fnum!=-1) { if (fnum!=-1) {
/* cout << " **Incomplete frame number " << fnum << " pnum " << pnum << " " << getFrameNumber(p) << endl;*/ /* cout << " **Incomplete frame number " << fnum << " pnum " << pnum << " " << getFrameNumber(p) << endl;*/