diff --git a/AnalysisADNE.C b/AnalysisADNE.C index b16a40d94ca473e6c314e71954aa8bc826cd0b0c..19436ffab4a6e533746590b1e80ef0ba41121909 100644 --- a/AnalysisADNE.C +++ b/AnalysisADNE.C @@ -77,7 +77,7 @@ int main(){ cerr<< "File Open"<<endl; printf("\033[32mInfo:: File Open\033[m \n"); //schaine.ReplaceAll("/space/clement/space/AGATA/OutMerger/",""); - //schaine.ReplaceAll("/data/exogamX/e826/acquisition/run/",""); + // schaine.ReplaceAll("/data/exogamX/e826/acquisition/run/",""); schaine.ReplaceAll("/data/e826X/e826/acquisition/run/",""); //schaine.ReplaceAll("/data/lise_zddX/test/acquisition/run/",""); if(schaine.BeginsWith("run_")){ diff --git a/GUser.C b/GUser.C index 7a57c72427787318cc681994e0a709ee123feffc..4d1f985efa3409d9342082ed9ca49d46b31c1691 100644 --- a/GUser.C +++ b/GUser.C @@ -101,7 +101,7 @@ GUser::GUser (GDevice* _fDevIn, GDevice* _fDevOut) fExogam2->ActivateGOCCETrack(false); //activate a second level of checking for GOCCE - fExogam2->SetPromptGate(0,64000); + fExogam2->SetPromptGate(0,1e5); fExogam2->SetCloverPosition(0,0,141.); //ecc# flange# dis[mm] // from 2022 and numexo2 eccId==flangeId @@ -192,8 +192,8 @@ void GUser::InitUser() if(fExogam2->IsCloverActive(k)){ MySpectraList->AddSpectrum(fExogam2->fMyHistoCloverECal[k],"Exogam/Clover_ECal"); MySpectraList->AddSpectrum(fExogam2->fMyHistoCloverECalACAdd[k],"Exogam/Clover_ECalACAdd"); - MySpectraList->AddSpectrum(fExogam2->fMyHistoCloverECalACAdd_TC[k],"Exogam/Clover_ECalACAdd_TC"); - MySpectraList->AddSpectrum(fExogam2->fMyHistoCloverECalACAdd_TQ[k],"Exogam/Clover_ECalACAddDC_TQ"); + MySpectraList->AddSpectrum(fExogam2->fMyHistoCloverECalACAdd_TC[k],"Exogam/Clover_ECalACAdd_TC"); //AddBack + AC + prompt /clover + MySpectraList->AddSpectrum(fExogam2->fMyHistoCloverECalACAdd_TQ[k],"Exogam/Clover_ECalACAddDC_TQ");//AddBack + AC /clover MySpectraList->AddSpectrum(fExogam2->fMyHistoCloverECalACAddRejectF[k],"Exogam/Clover_ECalACAddRejF"); MySpectraList->AddSpectrum(fExogam2->fMyHistoCloverTCal[k],"Exogam/Clover_TCal"); @@ -227,8 +227,14 @@ void GUser::InitUser() MySpectraList->AddSpectrum(fExogam2->fMyHistoGOCCENet,"Exogam/PSA") ; MySpectraList->AddSpectrum(fExogam2->fMyHistoGOCCEMirror,"Exogam/PSA") ; + MySpectraList->AddSpectrum(fExogam2->fMyHistoRMirror,"Exogam/PSA") ; MySpectraList->AddSpectrum(fExogam2->fMyHistoPhiMirror,"Exogam/PSA") ; - MySpectraList->AddSpectrum(fExogam2->fMyHistoSumECCE,"Exogam/Sum"); + MySpectraList->AddSpectrum(fExogam2->fMyHistoPSASurface,"Exogam/PSA") ; + MySpectraList->AddSpectrum(fExogam2->fMyHistoPSASurfaceCarte,"Exogam/PSA") ; + MySpectraList->AddSpectrum(fExogam2->ShortTrace,"Exogam/PSA") ; + MySpectraList->AddSpectrum(fExogam2->fMyHistoPSAChi2,"Exogam/PSA") ; + MySpectraList->AddSpectrum(fExogam2->fMyHistoPSAChi2_Radius,"Exogam/PSA") ; + MySpectraList->AddSpectrum(fExogam2->fMyHistoSumECCE,"Exogam/Sum"); //sumCore MySpectraList->AddSpectrum(fExogam2->fMyHistoSumECCT,"Exogam/Sum"); MySpectraList->AddSpectrum(fExogam2->fMyHistoSumGOCCEE,"Exogam/Sum"); @@ -239,7 +245,7 @@ void GUser::InitUser() MySpectraList->AddSpectrum(fExogam2->fMyHistoPatternGOCCEE,"Exogam/Exogam_Pattern"); MySpectraList->AddSpectrum(fExogam2->fMyHistoPatternGOCCET,"Exogam/Exogam_Pattern"); MySpectraList->AddSpectrum(fExogam2->fMyHistoECCcopy_val,"Exogam"); - MySpectraList->AddSpectrum(fExogam2->fMyHistoSumExogam2,"Exogam/Sum"); + MySpectraList->AddSpectrum(fExogam2->fMyHistoSumExogam2,"Exogam/Sum"); //AddBack + AC rejected MySpectraList->AddSpectrum(fExogam2->fMyHistoSumExogam2DC,"Exogam/Sum"); MySpectraList->AddSpectrum(fExogam2->fMyHistoSumExogam2FoldCond,"Exogam/Sum"); MySpectraList->AddSpectrum(fExogam2->fMyHistoSumExogam2DCFoldCond,"Exogam/Sum"); diff --git a/ListRun.txt b/ListRun.txt index 208d56ea9f94f2fda582caf56f58be466715832c..280f0fbbc1ee6012e7978154445b78102e109664 100644 --- a/ListRun.txt +++ b/ListRun.txt @@ -11,4 +11,5 @@ #exogamwithmirror #/data/exogamX/e826/acquisition/run/run_0041.dat.08-02-22_14h17m46s #exogamparisfull -/data/e826X/e826/acquisition/run/run_0072.dat.25-02-22_12h07m58s +#/data/e826X/e826/acquisition/run/run_0072.dat.25-02-22_12h07m58s +/data/e826X/e826/acquisition/run/run_0075.dat.25-02-22_17h20m44s diff --git a/TExogam2.cxx b/TExogam2.cxx index f4c1a26978e90e8c4af8abd09719ca593e6806db..3f12db898cc539a1562ded7b4ef0fc1254fdce6a 100644 --- a/TExogam2.cxx +++ b/TExogam2.cxx @@ -16,6 +16,12 @@ TExogam2::TExogam2(bool bspec) Goccetrack=GOCCEActive=ESSActive=false; LUTBool=false; PrevTS=0; + + + trace = new TF1("trace","[0]/(1.+exp((x-[1])/[2]))",-1000,1000); //wood-saxon like traces + traceExp = new TGraph(6); + + } @@ -167,7 +173,7 @@ if(BoolSpec){ HListExogam2.Add(fMyHistoECCECal[k][l-1]); sprintf(title,"ECC%d_%c_TCal",k,cristal[l-1]); - fMyHistoECCTCal[k][l-1]= new TH1F(title,title,4000,0,16000); + fMyHistoECCTCal[k][l-1]= new TH1F(title,title,4000,0,80000); //HListExogam2.Add(fMyHistoECCTCal[k][l-1]); sprintf(title,"ECC%d_%c_T30",k,cristal[l-1]); @@ -228,24 +234,34 @@ if(BoolSpec){ } } - fMyHistoSumECCE= new TH1F("SumEnergyCoreExogam2","SumEnergyCoreExogam2",4000,0,2000); + fMyHistoSumECCE= new TH1F("fMyHistoSumECCE","fMyHistoSumECCE",4000,0,2000); HListExogam2.Add(fMyHistoSumECCE); - fMyHistoSumECCT= new TH1F("SumTimeCoreExogam2","SumTimeCoreExogam2",4000,0,2000); + fMyHistoSumECCT= new TH1F("fMyHistoSumECCT","fMyHistoSumECCT",4000,0,2000); HListExogam2.Add(fMyHistoSumECCT); - fMyHistoSumGOCCEE= new TH1F("SumEnergySegExogam2","SumEnergySegExogam2",70000,0,70000); + fMyHistoSumGOCCEE= new TH1F("fMyHistoSumGOCCEE","fMyHistoSumGOCCEE",70000,0,70000); HListExogam2.Add(fMyHistoSumGOCCEE); - fMyHistoSumGOCCET= new TH1F("SumTimeSegExogam2","SumTimeSegExogam2",4000,0,2000); + fMyHistoSumGOCCET= new TH1F("fMyHistoSumGOCCET","fMyHistoSumGOCCET",4000,0,2000); HListExogam2.Add(fMyHistoSumGOCCET); fMyHistoGOCCENet= new TH1F("fMyHistoGOCCENet","fMyHistoGOCCENet",70000,0,70000); fMyHistoGOCCEMirror = new TH1F("fMyHistoGOCCEMirror","fMyHistoGOCCEMirror",70000,-35000,35000); - fMyHistoPhiMirror = new TH1F("fMyHistoPhiMirror","fMyHistoPhiMirror",1000,-200,200); + fMyHistoPhiMirror = new TH1F("fMyHistoPhiMirror","fMyHistoPhiMirror",400,0,400); + fMyHistoRMirror= new TH1F("fMyHistoRMirror","fMyHistoRMirror",100,0,100); + fMyHistoPSASurface= new TH2F("fMyHistoPSASurface","fMyHistoPSASurface",200,0,2000,1000,-200,200); + fMyHistoPSASurfaceCarte= new TH2F("fMyHistoPSASurfaceCarte","fMyHistoPSASurfaceCarte",200,-100,100,200,-100,100); + ShortTrace= new TGraph(6); + fMyHistoPSAChi2 = new TH1F("fMyHistoPSAChi2","fMyHistoPSAChi2",7000,0,700); + fMyHistoPSAChi2_Radius= new TH2F("fMyHistoPSAChi2_Radius","fMyHistoPSAChi2_Radius",7000,0,700,2000,0,1000); HListExogam2.Add(fMyHistoGOCCENet); HListExogam2.Add(fMyHistoGOCCEMirror); HListExogam2.Add(fMyHistoPhiMirror); + HListExogam2.Add(fMyHistoRMirror); + HListExogam2.Add(fMyHistoPSASurfaceCarte); + HListExogam2.Add(fMyHistoPSAChi2); + fMyHistoPatternECCE= new TH1F("PatternEnergyCoreExogam2","PatternEnergyCoreExogam2",200,0,200); HListExogam2.Add(fMyHistoPatternECCE); @@ -260,7 +276,7 @@ if(BoolSpec){ - fMyHistoSumExogam2= new TH1F("Exogam2AddB_AC_noDC","Exogam2AddB_AC_noDC",4000,0,2000); + fMyHistoSumExogam2= new TH1F("Exogam2AddB_AC_noDC","Exogam2AddB_AC_noDC",4000,0,2000);//AddBack + AC rejected HListExogam2.Add(fMyHistoSumExogam2); fMyHistoSumExogam2FoldCond= new TH1F("Exogam2AddB_AC_noDC_Fold1","Exogam2AddB_AC_noDC_Fold1",4000,0,4000); HListExogam2.Add(fMyHistoSumExogam2FoldCond); @@ -604,6 +620,9 @@ bool TExogam2::IsMFMExo(MFMExogamFrame *frame) double valf,valf2,valf3; bool result = false; bool debug = false; + //---------------PSA control + bool PSA_debug = false; + int EXO2BoardId, EXO2CrysId, MapFinger,MulNetCharge ; int statusSegment, binstatus; float s[4]={0.,0.,0.,0.}; @@ -613,11 +632,24 @@ bool TExogam2::IsMFMExo(MFMExogamFrame *frame) int q1,q2; float PSA_phi; float T30,T60,T90,PSA_r; - float I1,I2,I3; + TVector2 PSAHit(0,0); + float psa_T[6] = {0,0,0,0,0,0}; + float psa_Q[6] = {0,0,0,0,0,0}; + float psa_TD1[5] = {0,0,0,0,0}; + float psa_QD1[5] = {0,0,0,0,0}; + float psa_TD2[4] = {0,0,0,0}; + float psa_QD2[4] = {0,0,0,0}; + float QMAX, TMAX,chi2,Xp,Yp, Rp,ThetaP; + int Method=2; + float PSACutOFF=160.; + + //---------------PSA control //cout<<frame->GetTimeStamp()<<endl; clo=cri=MapFinger=-100; T30=T60=T90=PSA_r=PSA_phi=0; + QMAX= TMAX = -1e3; + if(LUTBool==false){ printf("\033[31m Error in TExogam2::IsMFMExo => The Look Up Table is not known: I won't be able to unpack the event \033[m \n"); } @@ -769,25 +801,24 @@ bool TExogam2::IsMFMExo(MFMExogamFrame *frame) MaxSeg[clo][cri]=seg; } }//end Loop Segement - /*if(MulNetCharge==1){ //PSA Studies - cout<<""<<endl; - cout<<"PSA Analysis || Single interaction only"<<endl; - printf("(%d)%f --- (%d)%f \n",psa[0],s[0],psa[1],s[1]); - printf("--- %f ---\n",valf); - printf("(%d)%f --- (%d)%f \n",psa[3],s[3],psa[2],s[2]); - cout<<"----------------------"<<endl; + if(MulNetCharge==1&&valf>100){ //PSA Studies + if(PSA_debug)cout<<""<<endl; + if(PSA_debug)cout<<"======================================="<<endl; + if(PSA_debug)cout<<"PSA Analysis || Single interaction only"<<endl; + if(PSA_debug)printf("(%d)%f --- (%d)%f \n",psa[0],s[0],psa[1],s[1]); + if(PSA_debug)printf("--- %f ---\n",valf); + if(PSA_debug)printf("(%d)%f --- (%d)%f \n",psa[3],s[3],psa[2],s[2]); + if(PSA_debug)cout<<"----------------------"<<endl; for(Int_t search = 0 ; search<4 ; search++){ if(psa[search]==0){ - PSANetcharge=search; + PSANetcharge=search; // search for the next charge of the event break; } } if(PSANetcharge==0){ q1=1; q2=3; - - } else if(PSANetcharge==3){ q1=0; @@ -799,19 +830,111 @@ bool TExogam2::IsMFMExo(MFMExogamFrame *frame) } PSA_phi=TMath::Log10(abs(s[q2])/abs(s[q1])); - cout<<"Phi " <<PSA_phi<<endl; - fMyHistoPhiMirror->Fill(0.26*PSA_phi*TMath::RadToDeg()+7.3); + PSA_phi=45.+(0.26*PSA_phi*TMath::RadToDeg()+10.); // in Deg 7.3 + + + + if(PSA_debug)cout<<"Phi " <<PSA_phi<<endl; + + psa_T[0]=Cal(T30,0.,-10.,0); + psa_T[1]= 0.; + psa_T[2]=-1.0*psa_T[0]; + psa_T[3]=Cal(T60,0.,10.,0); + psa_T[4]=Cal(T90,0.,10.,0); + psa_T[5]=psa_T[4]+400.; + + psa_Q[0]=0; + psa_Q[1]=0.; + psa_Q[2]=30.; + psa_Q[3]=60.; + psa_Q[4]=90.; + psa_Q[5]=100.; + if(PSA_debug)cout<<"Input Trace"<<endl; + for(Int_t psaStep = 0 ; psaStep<6; psaStep++){ + if(PSA_debug)cout<<psa_T[psaStep]<<" "<<psa_Q[psaStep]<<endl; + traceExp->SetPoint(psaStep,psa_T[psaStep],psa_Q[psaStep]); + } + + + if(Method==1){ + for(Int_t psaStep = 0 ; psaStep<5; psaStep++){psa_TD1[psaStep]=psa_T[psaStep+1];} + for(Int_t psaStep = 0 ; psaStep<5; psaStep++){psa_QD1[psaStep]=(psa_Q[psaStep+1]-psa_Q[psaStep])/(psa_T[psaStep+1]-psa_T[psaStep]);} + + for(Int_t psaStep = 0 ; psaStep<4; psaStep++){psa_TD2[psaStep]=psa_TD1[psaStep+1];} + for(Int_t psaStep = 0 ; psaStep<4; psaStep++){psa_QD2[psaStep]=(psa_QD1[psaStep+1]-psa_QD1[psaStep])/(psa_TD1[psaStep+1]-psa_TD1[psaStep]);} + + + for(Int_t psaStep = 0 ; psaStep<4; psaStep++){ + if(PSA_debug)cout<<" Time " << psa_TD2[psaStep] << " Charge "<< psa_QD2[psaStep]<<endl; + if(psa_QD2[psaStep]>QMAX){ + QMAX=psa_QD2[psaStep]; + TMAX=psa_TD2[psaStep]; + } + } + chi2=1; + + } + else if (Method==2){ + trace->SetParameter(0,100); + trace->SetParameter(1,100) ; + trace->SetParameter(2,-300); + for(Int_t trial = 0 ; trial < 5; trial++){ + traceExp->Fit("trace","Q"); + } + if(PSA_debug)cout<<"Chi2 is "<<trace->GetChisquare()<<endl; + + chi2=trace->GetChisquare(); + for(Int_t trial = psa_T[0] ; trial <psa_T[5] ; trial++){ + if(trace->Derivative2(trial)>QMAX){ + QMAX=trace->Derivative2(trial); + TMAX=1.*trial; + } + } + + + } + if(PSA_debug)cout<<"Max is "<<QMAX<<" at "<<TMAX<<endl; + + if(chi2<PSACutOFF){ //these are good traces + PSA_r=TMAX*0.05+0.0; + PSAHit.SetMagPhi(PSA_r,PSA_phi*TMath::DegToRad()); + //cout<<PSAHit.Phi()*TMath::RadToDeg()<<endl; + //cout<<"then rotate by "<<1.*PSANetcharge*(90.*TMath::DegToRad())<<endl; + + //PSAHit.Rotate((double)(PSANetcharge*90.*TMath::RadToDeg())); + Xp=PSAHit.X()*TMath::Cos(PSANetcharge*90.*TMath::RadToDeg())-PSAHit.Y()*TMath::Sin(PSANetcharge*90.*TMath::RadToDeg()); + Yp=PSAHit.X()*TMath::Sin(PSANetcharge*90.*TMath::RadToDeg())+PSAHit.Y()*TMath::Cos(PSANetcharge*90.*TMath::RadToDeg()); + Rp=sqrt(Xp*Xp+Yp*Yp); + + if(Xp>0&&Yp>=0)ThetaP=TMath::ATan(Yp/Xp); + if(Xp>0&&Yp<0)ThetaP=TMath::ATan(Yp/Xp)+2*TMath::Pi(); + if(Xp<0)ThetaP=TMath::ATan(Yp/Xp)+TMath::Pi(); + if(Xp==0&&Yp>0)ThetaP=TMath::Pi()/2.; + if(Xp==0&&Yp<0)ThetaP=3*TMath::Pi()/2.; + PSAHit.SetMagPhi(Rp,ThetaP); + + //cout<<PSAHit.Phi()*TMath::RadToDeg()<<endl; + if(BoolSpec)fMyHistoPhiMirror->Fill(PSAHit.Phi()*TMath::RadToDeg());// in Deg + + if(BoolSpec)fMyHistoRMirror->Fill(PSA_r); + if(BoolSpec)fMyHistoPSASurface->Fill(PSA_r,PSAHit.Phi()*TMath::RadToDeg()); + if(BoolSpec)fMyHistoPSAChi2_Radius->Fill(chi2,PSA_r); + + if(BoolSpec)fMyHistoPSAChi2->Fill(chi2); + if(BoolSpec)fMyHistoPSASurfaceCarte->Fill(PSAHit.X(),PSAHit.Y()); + - I1=(60.-30.)/(1.*T60-T30); - I2=(90.-60.)/(1.*T90-T60); - I3=(I2-I1)/(1.*T90-T60); - cout<<"I1 " <<I1<<" I2 "<<I2 <<endl; - }*/ + for(Int_t psaStep = 0 ; psaStep<6; psaStep++){ + if(PSA_debug)cout<<psaStep<<"---- psa_T / psa_Q "<<psa_T[psaStep]<< "/" <<psa_Q[psaStep]<<endl; + ShortTrace->SetPoint(psaStep,psa_T[psaStep],psa_Q[psaStep]); + } + } + } //------------------Now Anti Comption UShort_t sum = frame->ExoGetBGO()+frame->ExoGetCsi(); //do a sum of CsI+BGO energies - if(sum>10){ + if(sum>150){ NoComptonCore[clo][cri]=false; NoCompton[clo]=false; } @@ -841,6 +964,7 @@ float Theta_Gamma, Phi_Gamma; bool IsPrompt[16][4]; float SumCalorimeter; + SumCalorimeter=0; //------------------Time Treat for(Int_t c=0;c<16;c++){ diff --git a/TExogam2.h b/TExogam2.h index f3695a8bcfdf67c093f9071a74d397f617e1c955..c089bd531c3b51ea8048d289039c24fe4587b41d 100644 --- a/TExogam2.h +++ b/TExogam2.h @@ -14,11 +14,13 @@ #include <MFMExogamFrame.h> #include "TVector3.h" +#include "TVector2.h" #include "TF1.h" #include "TTree.h" #include "TH1.h" #include "TH2.h" #include "TMath.h" +#include "TGraph.h" #include <TObject.h> #include "TRandom.h" #include "TString.h" @@ -76,7 +78,17 @@ class TExogam2 : public TDetector { TH1F * fMyHistoGOCCENet ; TH1F * fMyHistoGOCCEMirror ; TH1F * fMyHistoPhiMirror ; - + TH1F * fMyHistoRMirror ; + TH2F * fMyHistoPSASurface ; + TH2F * fMyHistoPSASurfaceCarte ; + TGraph *ShortTrace; + TF1 *trace; + TGraph *traceExp; + + + TH1F * fMyHistoPSAChi2 ; + TH2F * fMyHistoPSAChi2_Radius ; + TH1F * fMyHistoESS_BGO[16][4]; TH1F * fMyHistoESS_CSI[16][4]; diff --git a/Utils/files.h b/Utils/files.h new file mode 100755 index 0000000000000000000000000000000000000000..f8fa62bfdab8d87bfd59b474fbc602117e7c680e --- /dev/null +++ b/Utils/files.h @@ -0,0 +1,266 @@ +#ifndef FILES_H +#define FILES_H 1 + +#define matsize 4096 +#define io_read 0 +#define io_write 1 + + +class files{ + int file_id; + char filename[80]; + char text[100]; + int lum; + FILE *fp; + int fd; +public: + char EUROGAM_NAME[80]; + files(); + files(const char*); + struct ganil_header{ + unsigned long num_run; // ! numero run (a la sauvegarde) + char reserve_1[12]; // ! + char nom_run[16]; // ! nom run (facultatif, 15 car utiles) + long int num_spectre; // ! numero spectre + char type_spectre[2]; // ! type ("1D" ou "2D") + char reserve_2[10]; // ! + char nom_spectre[16]; // ! nom du spectre (15 car utiles) + char date[12]; // ! date de creation ("jj-mmm-aaaa") + char heure[8]; // ! heure de creation ("hh:mn:ss") + long int nat_spectre; // ! reserve ulterieure (=0) + unsigned long codeur_x; // ! val max du codeur en x + unsigned long codeur_y; // ! val max du codeur en y + unsigned long taille_canal; // ! taille utile par canal (16 ou 32) + char type_canal[4]; // ! type variable ("I*4" , "I*2") + unsigned long dim_x; // ! dimension du spectre en x + unsigned long dim_y; // ! dimension du spectre en y + unsigned long min_x; // ! numero du canal min en x + unsigned long min_y; // ! numero du canal min en y + unsigned long max_x; // ! numero du canal max en x + unsigned long max_y; // ! numero du canal max en y + char commentaire[80]; // ! commentaire utilisateur + char nom_par_x[16]; // ! nom du parametre x (15 car utiles) + char nom_par_y[16]; // ! nom du parametre y (15 car utiles) + char unite_x[8]; // ! unite en x (8 car utiles) + char unite_y[8]; // ! unite en y (8 car utiles) + }header_ganil; + + struct header_ans{ + char revision_code[2]; //0-1 + char file_name[12]; //2-13 + char spectrum_id[72]; //14-85 + char Start_Date[8]; //86-93 + char Stop_Date[8]; //94-101 + char ElapsedRT[8]; //102-109 + char ElapsedLT[8]; //110-117 + char ElapsedPass[4]; //118-121 + char DeadTime[4]; //122-125 + char Hard_Setup_AMP[42]; //126-167 + char Hard_Setup_ADC[84]; //168-251 + char Hard_Setup_Bias[30]; //252-281 + char Ener_Calib_Data[32]; //282-313 + char Shape_Calib_Data[24]; //314-337 + char Eff_Calib_Data[48]; //338-385 + char Detector[208]; //386-593 + char Hardware_spec[24]; //594-617 + char Hardware_setup_Presets[20]; // 618-637 + char Tool_Settings[28]; //638-665 + char Background[24]; // 666-689 + char Microanalysis[48]; // 690-737 + char Oxford[10]; // 738-747 + char Sample[20]; // 748-767 + char More_Tools[20]; // 768-787 + char Reserved[226]; // 788-1013 + char SpecStart[2]; + char SpecEnd[2]; + char NROIs[2]; + char AcqSysType[2]; + char AcqSysLen[2]; + }header_ans; + + struct in2p3_header{ + char nomrun[16]; + int numrun; + char nomspe[16]; + int nuspe; + char date[20]; + int type; + int format; + int nbcanx; + int candebx; + int nbcany; + int candeby; + int nbcanz; + int candebz; + int nbcant; + int candebt; + char oct[4]; + char idaplic[8]; + int hlong; + char applic[16]; + char version[8]; + }header; + + struct header_eurogam + { + unsigned int magic_number; + unsigned int version; // expect 1 !!! + char spectrum_name[32]; // spectrum name + unsigned int fold; + char creation_date[20]; //creation date + char modification_date[20]; // modif date + + unsigned int base_info[8]; // base info for dim 1-8 + + unsigned int range[8]; // range info for dim 1-8 + + unsigned int fip[32]; // information pointer + int anotation_pointer[8]; + int calibration_pointer[8]; + int efficiency_pointer[8]; + unsigned int data_array_descriptor1[5]; + // 32 bits : array layout (0 for histo/matrix) + // 32 bits : array type (5 for 32 bit signed ) + // 32 bits : reserved for future + // 32 bits : reserved for future + // 32 bits : pointer + unsigned int data_array_descriptor2[5]; + int base_address_string_space; + int string_free_space; + unsigned int top_string_space; + unsigned int base_address_count_space; + unsigned int count_free_space; + unsigned int top_count_space; + unsigned int unused[20]; // not used ... + }header_eurogam; + + void check_little_endian(void); + unsigned long swap(unsigned long); + void swap4b(char *); + void swap2b(char *); + int open_file(void); + int open_file_r(void); + int open_file_r(char *); + int fopen_file(void); + void get_name(const char *); + int is_file_here(void); + int is_file_here(char *); + void remove(void); + void remove(char *); + void touch(void); + void put_ex_mat(); + void put_ex_mat(char *); + void put_ex_spe(); + void put_ex_spe(char *); + void put_ex_txt(char *); + void put_ex_ans(char *); + void put_ex_asc(char *); + void put_ex_s(); + void put_ex_s(char *); + void put_ex_mtr(); + int read_line(FILE *, char *); + void specdump(short *, int, FILE*, FILE*); + void specdump(int *, int, FILE*, FILE*); + void specdump(unsigned int *, int, FILE*, FILE*); + void specdump(float *, int, FILE*, FILE*); + void specdump(int **, int, int, FILE*, FILE *); + void specdump(short **, int, int, FILE*, FILE *); +// for radware .mat matrices + int write_mat_rad(short *, FILE *, FILE *); + int write_mat_rad(short *, char *, FILE *, FILE *); + void write_mat_rad(int *); + int read_mat_rad(short *, FILE *, FILE *); + int read_mat_rad(short *, char *, FILE *, FILE *); +// for mktri matrices format + int read_mktri(short *, int); + int write_mktri(short *, int); + int check_mat(int, int*, FILE*, FILE* ); + int check_mat(int, FILE*, FILE* ); +// for in2p3 format + int read_header_in2p3(); + void write_header_in2p3(FILE *); + int read_sp_in2p3(unsigned int *); + int read_sp_in2p3(unsigned int **); + int get_type_in2p3(){return(header.type);}; + int get_nx_in2p3(){return(header.nbcanx);}; + int get_ny_in2p3(){return(header.nbcany);}; + char* get_nom_in2p3(){return(header.nomspe);}; +// for gpsi format + int get_nsp_gpsi(char *); + int get_nx_gpsi(int, char *); + int get_ny_gpsi(int, char *); + void myrgspe(char *, int, float *,int , int, char *, int); + void myrgspe(FILE *,char *, int, float *,int , int, char *, int); +// for cambda format + int get_size_cambda(char *); + void get_cambda(char *, float *, int); + void get_att_mtr(void); + int read_mtr(float **); +// for ascii format + int get_size_ascii(char *); + void get_ascii(char *, float *, int); + void write_tab(float *, int); + void write_tab(int *, int); +// for eurogam format + int swaping_eurogam; + int open_eurogam(char *); + void close_eurogam(void); + int get_parameters_eurogam(char *, int *, int *, int *, FILE *, FILE *); + int get_header_eurogam(char *, FILE *, FILE *); + int get_header_eurogam(char *); + void print_header_eurogam(void); + int check_magicnumber_eurogam(char *, FILE *, FILE *); + int check_magicnumber_eurogam(char *); + void swap_header_eurogam(void); + + int get_fold_eurogam(void){return header_eurogam.fold;}; + int get_size_eurogam(void){return header_eurogam.count_free_space;}; + int get_size1_eurogam(void){return header_eurogam.range[0];}; + int get_size2_eurogam(void){return header_eurogam.range[1];}; + int get_type_eurogam(void){return + header_eurogam.data_array_descriptor1[1];}; + int get_tab1D_eurogam(unsigned char *); + int get_tab1D_eurogam(char *); + int get_tab1D_eurogam(unsigned short *); + int get_tab1D_eurogam(short *); + int get_tab1D_eurogam(unsigned int *); + int get_tab1D_eurogam(float *); + int get_tab1D_eurogam(int *); +// int get_tab2D_eurogam(short *); + int get_tab2D_eurogam(short *, int, int); + int get_tab2D_eurogam(short **, int, int); +// int get_tab2D_eurogam(int *, int); + int get_tab2D_eurogam(int *, int, int); + int get_tab2D_eurogam(int **, int, int); + int init_header_eg(int, int , const char *, int, FILE *, FILE *); + int write_eg(int*, int, int , const char *, FILE *, FILE *); + int write_eg(int*, int, int , const char *); + int write_eg(int**, int, int , const char *, FILE *, FILE *); + int write_eg(unsigned int*, int, int , const char *, FILE *, FILE *); + int write_eg(short** , int, int , const char *, FILE *, FILE *); + int write_eg(float*, int, int , const char *, FILE *, FILE *); +// for ganil format + int get_header_ganil(char *); + void print_header_ganil(void); + int get_dimx_ganil(void){return swap(header_ganil.dim_x);}; + int get_ganil_integer4(unsigned long *); +// for radware monodim format + int write_spe(char *, float *, int, FILE *, FILE *); + int read_spe(float *, char *, int *); + int read_spe(char *,float *, char *, int *, FILE *, FILE *); +// for xmgr + int tab_to_xmgr(int, float*, float*, float*); + void write_xmgr(int, float *, float *); +// for ans format + int get_header_ans(char *); + int open_ans(char *); + void close_ans(void); + void print_header_ans(void); + int get_size_ans(); + int read_sp_ans(float *); +// for asc format + int read_header_asc(char *); + void read_asc(float *, char *, int); +}; + +#endif diff --git a/Utils/psa.txt b/Utils/psa.txt new file mode 100644 index 0000000000000000000000000000000000000000..58643781913c1962b06741c8d3a5da7029187750 --- /dev/null +++ b/Utils/psa.txt @@ -0,0 +1,6 @@ +-147.197 0 +0 0 +147.197 30 +279.857 60 +484.151 90 +884.151 100 diff --git a/Utils/root2spe b/Utils/root2spe new file mode 100755 index 0000000000000000000000000000000000000000..1f7b2605c3276f66172204faeecadd5a36693c70 Binary files /dev/null and b/Utils/root2spe differ diff --git a/Utils/root2spe.C b/Utils/root2spe.C new file mode 100755 index 0000000000000000000000000000000000000000..18fed1ea4bd6a4ee906034cc67fe42bbd36d1799 --- /dev/null +++ b/Utils/root2spe.C @@ -0,0 +1,264 @@ +#include <string.h> +#include <stdio.h> +#include <iostream> +#include <fstream> +#include <stdlib.h> +#include <math.h> +#include <time.h> +#include <ctype.h> +#include "TStyle.h" +#include "TCanvas.h" +#include "TGraph.h" +#include "TH2F.h" +#include "TChain.h" +#include "TTree.h" +#include "TNtuple.h" +#include "TFile.h" +#include "TRandom.h" +#include "TVector3.h" +#include "TRotation.h" +#include "TLorentzVector.h" +#include "TGenPhaseSpace.h" +#include "TMath.h" +#include "TRint.h" +#include "TString.h" +#include "TKey.h" +#include "files.h" + +/* root2spe spectrum converter root->spe + 1D histograms only + Ch. Theisen + First Version : Feb 2005 + + compile with : +g++ root2spe.C -o root2spe `$ROOTSYS/bin/root-config --libs` `$ROOTSYS/bin/root-config --cflags` +*/ + + +class files File; + +files::files(){ + + strcpy(EUROGAM_NAME,"Eurogam spectrum"); + +} + +files::files(const char *fname){ + + strcpy(EUROGAM_NAME,fname); + +} +int files::is_file_here(char *file){ +FILE *fpa; + + fpa = fopen(file,"r"); + if (fpa == NULL) return(0); + fclose(fpa); + return(1); +} + +void files::swap4b(char *buf){ +char c; + + c = buf[3]; buf[3] = buf[0]; buf[0] = c; + c = buf[2]; buf[2] = buf[1]; buf[1] = c; + +} + +int files::write_spe(char *file, float *spectre, int size, // output name , spectra array[nbin], nbin, + FILE *stream1, FILE *stream2){ +FILE *fprad; +int dum; +int j = 24; +int toto; +int tab[4] = {0,1,1,1}; +int i; +float *newspec; +int newsize; + + + newsize = size; + if (size > 16384){ + newsize=16384; + printf("Initial size = %d; resizing to %d \n",size,newsize); + } + + tab[0] = newsize; + toto = newsize * sizeof(float); + + fprad = fopen(file,"w"); + if (fprad == NULL){ + fprintf(stream1,"Cannot write %s \n",file); + if (stream2) + fprintf(stream2,"Cannot write %s \n",file); + fclose(fprad); + return(0); + } + if (LITTLE_ENDIAN == 1){ + swap4b((char*)&j); + for(i=0; i<4; i++) + swap4b((char*)&(tab[i])); + swap4b((char*)&toto); + newspec = new float[newsize]; + for(i=0; i<newsize; i++){ + + newspec[i] = spectre[i]; + swap4b((char*)&(newspec[i])); + } + } + + dum = fwrite(&j, sizeof(int), 1, fprad); + dum = fwrite(filename, 8, 1, fprad); + dum = fwrite(tab, sizeof(int), 4, fprad); + dum = fwrite(&j, sizeof(int), 1, fprad); + + dum = fwrite(&toto, sizeof(int), 1, fprad); + if (LITTLE_ENDIAN == 1){ + dum = fwrite(newspec, sizeof(float), newsize, fprad); + delete [] newspec; + } + else + dum = fwrite(spectre, sizeof(float), newsize, fprad); + + dum = fwrite(&toto,sizeof(int),1,fprad); + if (fprad == NULL){ + fprintf(stream1,"Cannot write %s \n",file); + if (stream2) + fprintf(stream2,"Cannot write %s \n",file); + + fclose(fprad); + return(0); + } + fclose(fprad); + + return(1); +} + + +void files::put_ex_spe(char *name){ +int len; + len = strlen(name); + if ( ( name[len-1] == 'e' ) && + ( name[len-2] == 'p' ) && + ( name[len-3] == 's' ) && + ( name[len-4] == '.' ) ) return; + sprintf(name,"%s.spe",name); +} +int process(TH1F *histo, int option){ +char name[80]; +char radname[80]; +char type[80]; +int nx,ny; +int binning,range,start; +float *spec; +float ch,st; +int i,j; + + sprintf(name,histo->GetName()); + printf("processing %s ...\n",name); + sprintf(type,histo->ClassName()); + if (strncmp(type,"TH1",3) != 0){ + printf(" %s is not a TH1 \n",name); + return(1); + } + + nx = histo->GetNbinsX(); + ny = histo->GetNbinsY(); + printf(" Nombre de bins X : %d \n",nx); + printf(" Nombre de bins Y : %d \n",ny); + printf(" First bin : %f \n",histo->GetBinLowEdge(1)); + printf(" Last bin : %f \n",histo->GetBinLowEdge(nx)); + + spec = new float[nx]; + for(i=1; i<=nx; i++){ + ch = histo->GetBinLowEdge(i); + st = histo->GetBinContent(i); + //printf("%d %f %f \n",i,ch,st); + spec[i-1] = histo->GetBinContent(i); + } + + + + switch(option){ + case 1: + sprintf(radname,name); + File.put_ex_spe(radname); + printf(" writing %s \n",radname); + break; + case 2: + printf("Enter output .spe filename : "); + scanf("%s",radname); + File.put_ex_spe(radname); + break; + } + + File.write_spe(radname, spec, nx, stdout, NULL); // output name , spectra array[nbin], nbin, + delete [] spec; + return(1); +} + +int main(int argc, char *argv[]){ +int mode; // 1 for only one 2 for all +int dum; +char file[80]; +char specname[80]; +TH1F *histo; + + printf("--------------------------------------------------- \n"); + printf(" root2spe : a program to convert root spectrum \n"); + printf(" to Radford spectrum format \n"); + printf(" use root2spe -all to convert all histos \n"); + printf("--------------------------------------------------- \n"); + + if (argc == 1){ + mode = 1; + } + else{ + dum = strcmp(argv[1],"-all"); + if (dum != 0){ + printf("unvalid option %s \n",argv[1]); + exit(0); + } + mode = 2; + } + + printf("root file : "); + scanf("%s",file); + + if (!File.is_file_here(file)){ + printf("This file does not exists ...\n"); + exit(0); + } + + + TFile *file1 = new TFile(file); + if (file1->IsZombie()){ + printf("%s is not a root file ... \n",file); + exit(0); + } + + switch(mode){ + case 1: + file1->ls(); + printf("histogram name : "); + scanf("%s",specname); + + histo = (TH1F*)file1->Get(specname); + if (histo == NULL){ + printf("Error : %s does not exist... \n",specname); + exit(0); + } + process(histo,2); + break; + case 2: + TIter nextkey(file1->GetListOfKeys()); + TKey *key; + while(key = (TKey*)nextkey()){ + TH1F *histo = (TH1F*)key->ReadObj(); + process(histo,1); + } + break; + } + file1->Close(); + +}