Commit 90763052 authored by Jérémie Dudouet's avatar Jérémie Dudouet
Browse files

Merge branch 'preprod' into 'master'

Integration of FIPPSSpy and New version of Cubix, including radware's tools

See merge request IPNL_GAMMA/gammaware!45
parents c6ba9341 09e3e3de
Pipeline #39132 passed with stage
in 7 minutes and 15 seconds
......@@ -57,6 +57,21 @@ link_directories(${CMAKE_LIBRARY_OUTPUT_DIRECTORY})
#
# configure file
#
#TEST current configuration OS
if ( UNIX )
set(CMAKE_OS_TYPE 1)
endif()
if ( WIN32 )
set(CMAKE_OS_TYPE 2)
endif()
if ( CYGWIN )
set(CMAKE_OS_TYPE 3)
endif()
if ( APPLE )
set(CMAKE_OS_TYPE 4)
endif()
configure_file (
"${CMAKE_CURRENT_SOURCE_DIR}/GwConfig.h.cmake"
"${CMAKE_BINARY_DIR}/include/GwConfig.h"
......
......@@ -43,3 +43,14 @@
#define GW_MacrosPath "@GW_DATA_DIR@/macros";
#endif
//////////////////////////////
////// Define OS_TYPE ////////
//////////////////////////////
#define OS_LINUX 1
#define OS_WINDOWS 2
#define OS_CYGWIN 3 // rather useless ==> to be removed ??
#define OS_APPLE 4
#define OS_TYPE @CMAKE_OS_TYPE@
......@@ -171,6 +171,25 @@ BaseGEM::BaseGEM():
fAD.SetOwner(kTRUE);
}
BaseGEM::BaseGEM(LevelScheme *lev) :
TObject(),
fLevelScheme(lev), fFeeding(), fRandFeeding(), fRandDown(), fAD(), fRand0(), fTmp(), fEoC(),
fDirection(kDown),
fMinEnergyFactor(10.0),
fADType(0),
fCutLifeTime(0.0000000001), /* 0.1*ns from Krane 1971 */
fSigma(0.0)
{
fTmp.SetOwner(kTRUE);
fRandDown.SetOwner(kTRUE);
fFeeding.SetOwner(kTRUE);
fRandFeeding.SetOwner(kTRUE);
fAD.SetOwner(kTRUE);
InitRandom();
}
BaseGEM::~BaseGEM()
{
BaseGEM::Clear(""); // clear everything
......@@ -576,8 +595,8 @@ Int_t BaseGEM::DoCascade(TSeqCollection &cas, TSeqCollection &directions, Int_t
obj = ((RandObj *)fRandFeeding.At( which_first ))->Rand();
while ( obj != &fEoC ) {
link = (Link *)fLevelScheme->GetLinks().At(obj->GetUniqueID());
result_one_link += link->DoCascade(cas,o);
result_one_link = link->DoCascade(cas,o);
for (Int_t i = result; i < result + result_one_link; i++) { // make sur the direction collection is synchronize with cas
if ( directions.At(i) == 0 ) {
directions.AddAt(new TLorentzVector(),i);
......
......@@ -141,7 +141,8 @@ namespace Gw {
virtual void DoAngularDistribution(Int_t which_gamma, TLorentzVector &, Bool_t forceiso = false);
public:
BaseGEM();
BaseGEM();
BaseGEM(LevelScheme *);
virtual ~BaseGEM();
//! To get all the entry point in this level scheme with their intensities
......
......@@ -6,6 +6,7 @@
# Add modules
#
add_subdirectory (cubix)
add_subdirectory (fipps_spy)
......
......@@ -38,7 +38,7 @@ set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++11 -pthread -Wno-unused-paramete
# root
INCLUDE_DIRECTORIES( ${ROOT_INCLUDE_DIR} )
LINK_DIRECTORIES( ${ROOT_LIBRARY_DIR} )
SET(EXTRA_EXTERNAL_LIBRARIES ${EXTRA_EXTERNAL_LIBRARIES} Matrix)
SET(EXTRA_EXTERNAL_LIBRARIES ${EXTRA_EXTERNAL_LIBRARIES} Matrix Rint)
# gw
set(EXTRA_INTERNAL_LIBRARIES GWPHYSICS )
......
ENCAL OUTPUT FILE
TMEUHS.STO EU152.SOU TM 159 EU152 HIGH SPIN CALIBRATION
2 0.000000000000D+00 1.000000000000D+00 0.000000000000D+00
0.000000000000D+00 0.000000000000D+00 0.000000000000D+00
### Cube file
CubeFileName MyCube.cub
### 2D projection file
2DProj MyCube.2dp
### LUT file
LUT MyCube.tab
### Contraction factor (default 1)
CompressFact 1
### Total projection file
TotalProj MyCube.spe
### Backgroud file
BackgroudFile MyCubeBG.spe
### Calibration file (not mandatory)
#CalFile MyCube.aca
### Efficiency file (not mandatory)
#EffFile MyCube.aef
#######################
### Cube Properties ###
#######################
# Name of the 3D cube file
CubeFileName MyCube.cub
# name of the look-up table file to be used to map ADC chs to cube chs
TabFile test.tab
# Size to be used for the "scratch" disk file, in MB
ScratchSize 512
# Compression factor (to obtain binning <1 keV, one needs to scale the energy)
# For example, to obtain a minimal 0.5 KeV/Channel, we fill the Cube on 8k channel for 4keV range ==> Compression factor = 2
# Must be similar to the value in the GenBinning conf file
CompressionFactor 1
################################
## Input Root Tree Properties ##
################################
#Name of the input TTree (regular expressions can be used as for building TChain)
TreeFile MyTree.root
#Name of the TTree
TreeName TreeName
#Events to process (0 => All)
NEvents 1000
#Name of the branch corresponding to the gamma-ray energy
EGamma EGamma
#Name of the branch corresponding to the gamma-ray muliplicity
GammaMult GammaMult
#Cut to be apllied (see TTreeFormulaExpression)
Cut None
EFFIT PARAMETER FILE 19-FEB-92 10:14:04
0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00
0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+02 0.00000000E+03
#Name of the tab file
TabFile test.tab
#To fit the FWHM vs Energy function:
#Fit function: sqrt([0]*[0] + [1]*[1]*x/1000. + [2]*[2]*x/1000.*x/1000.)
#For each new point: FWHM Energy (keV) Value (keV)
FWHM 122. 1.35
FWHM 1333. 2.35
FWHM 2500. 3.00
# Compression factor (to obtain binning <1 keV, one needs to scale the energy)
# For example, to obtain a minimal 0.5 KeV/Channel, we fill the Cube on 8k channel for 4keV range ==> Compression factor = 2
CompressionFactor 1
#Number of bin per FWHM
Precision 1
#RangeMax (keV)
RangeMax 8192
ROOTCFLAGS := $(shell root-config --cflags)
ROOTLIBS := $(shell root-config --libs)
CXX = g++
CXXFLAGS = -O3 -I . -w
CXXFLAGS = -O3 -I ./src -w -g $(ROOTCFLAGS) -Wno-write-strings -Wno-unused-result
LDFLAGS = ${ROOTLIBS}
obj_CubeBuilder = src/CubeBuilder.o src/progressbar.o
obj_GenBinning = src/GenBinning.o src/util.o
ALL = CubeBuilder GenBinning
.PHONY: all
all: $(ALL)
CubeBuilder: $(obj_CubeBuilder)
$(CXX) $^ -o $@ $(LDFLAGS) $(CXXFLAGS)
GenBinning: $(obj_GenBinning)
$(CXX) $^ -o $@ $(LDFLAGS) $(CXXFLAGS)
.PHONY: clean
clean:
rm -f $(ALL) src/*.o
ReadMe CubeBuilder:
1 - Generate the adaptative binning file:
./GenBinning ConfBinning.txt
This will produce:
- The binning file: Binning.tab
- A width spectrum file: Binning.spe
2 - Generate the Cube
./CubeBuilder ConfBuildCube.conf
This will produce:
- The cube file: MyCube.cub
3 - Generate the total projection
(radware prog) pro3d => MyCube.cub
This will produce:
- The total 2D projection spectrum: MyCube.2dp
- The total 1D projection spectrum: MyCube.spe
4 - Generate the background (in gf3)
in gf3:
4-1: Open the 1D projection: MyCube.spe
4-2: Divide by the width spectrum: DV, then Binning.spe
4-3: Calculate background: BG, then ok, then click limits on the graph
4-4: Multiply by the width spectrum: MS, then Binning.spe
4-5: Save the background spectrum: WS, then file name
This diff is collapsed.
#define RW_MAXFILES 9999 /* Max files per tape to read */
#define RW_MAXRECORDS 1000000000 /* Max number of records per tape to read */
#define RW_DATFILE "incub8r.inc" /* Defaults file for hypercube */
#define RW_CHNK_S 32 /* size of update chunk, in units of 1/4 MB */
/*
RW_LB2 = ceil((float)(ceil((float)((RW_MAXCH/8)*((RW_MAXCH/8)+1)*((RW_MAXCH/8)+2)/6)/128))/RW_CHNK_S)
RW_LB1 = RW_LB2*RW_CHNK_S
*/
//for 1400 channels
//#define RW_MAXCH 1400 /* Maximum channels on cube axis */
// /* (must be multiple of 8) */
//#define RW_LB2 222 /* length required for buf2, depends on RW_LB1 */
// /* RW_LB2 >= RW_LB1/RW_CHNK_S */
//#define RW_LB1 7104 /* length required for buf1, depends on RW_MAXCH */
// /* RW_LB1 > x(x+1)(x+2)/(6*128), */
// /* where x = RW_MAXCH/8 */
// /* must also be a multiple of RW_CHNK_S */
//for 4096channels
//#define RW_MAXCH 4096
//#define RW_LB1 175808
//#define RW_LB2 5494
//for 8196 channels
#define RW_MAXCH 8192
#define RW_LB1 1402208
#define RW_LB2 43819
#define RW_DB1 90 /* depth of buf1 */
#define RW_DB2 180 /* depth of buf2 not to exceed RW_SCR_RECL */
/* byte records in scr file */
#define RW_SCR_RECL 32768 /* Record length of scr file */
/* must be >= RW_DB2*2*(RW_DB1+1) + 4 */
#define RW_MAXMULT 40
/************************************************
* 3D CUBE FILE FORMAT
* 1024byte file header
* 4096byte data records
*
* each data record contains:
* variable number of bit-compressed 8x8x4 mini-cubes
* 8bytes header, 4088 bytes of data
*/
/* 3d cube file header */
typedef struct {
char id[16]; /* "Incub8r3/Pro4d " */
int numch; /* number of channels on axis */
int bpc; /* bytes per channel, = 4 */
int cps; /* 1/cps symmetry compression, = 6 */
int numrecs; /* number of 4kB data records in the file */
char resv[992]; /* FUTURE flags */
} FHead3D;
/* 3Drecord header */
typedef struct {
int minmc; /* start minicube number, starts at 0 */
short nummc; /* number of minicubes stored in here */
short offset; /* offset in bytes to first full minicube */
} RHead3D;
typedef struct {
RHead3D h;
unsigned char d[4088]; /* the bit compressed data */
} Record3D; /* see the compression alg for details */
#include <stdio.h>
#include <stdlib.h>
#include <unistd.h>
#include <math.h>
#include <string.h>
#include <fstream>
#include "util.h"
#include "TGraph.h"
#include "TF1.h"
#include "TString.h"
#include "TCanvas.h"
#include "TApplication.h"
#include "TSystem.h"
using namespace std;
float wid_spec[16384];
float wa = 3.0f, wb = 1.0f, wc = 4.0f, cfactor = 1.0f, chs = 2.0f;
float w, x, f1, f2, f3, wa2, wb2, wc2;
int iclo = 200, ichi = 5000;
int i1, i2, i3, hi, lo, nc, i, j, ich2;
int nclook, lookmin, lookmax, rl = 0;
short looktab[16384];
char buf[32];
TString TabName = "FWHM.tab";
Bool_t fFixedBinning = false;
bool CalcFWHMParameters(vector< pair<float,float> > calibpoints) {
TApplication *theApp = new TApplication("App", 0, 0);
TGraph *CalibGraph = new TGraph;
CalibGraph->SetNameTitle("CalibFWHM","CalibFWHM");
TF1 *CalibFunc = new TF1("CalibFunc","sqrt([0]*[0] + [1]*[1]*x/1000. + [2]*[2]*x/1000.*x/1000.)");
for(uint i=0 ; i<calibpoints.size() ; i++) {
CalibGraph->SetPoint(CalibGraph->GetN(),calibpoints.at(i).first,calibpoints.at(i).second);
}
TCanvas *c = new TCanvas;
CalibGraph->SetMarkerStyle(20);
CalibGraph->SetMarkerSize(2);
CalibGraph->SetMarkerColor(kBlue);
CalibGraph->Draw("APL");
CalibGraph->Fit(CalibFunc,"SQ");
wa2 = CalibFunc->GetParameter(0)*CalibFunc->GetParameter(0);
wb2 = CalibFunc->GetParameter(1)*CalibFunc->GetParameter(1);
wc2 = CalibFunc->GetParameter(2)*CalibFunc->GetParameter(2);
cout<<" -- Fit Function: FWHM = sqrt(F*F + G*G*x + H*H*x*x) (x = ch.no. /1000)"<<endl;
cout<<" --> F = "<< CalibFunc->GetParameter(0)<<endl;
cout<<" --> G = "<< CalibFunc->GetParameter(1)<<endl;
cout<<" --> H = "<< CalibFunc->GetParameter(2)<<endl;
c->Update();
c->Modified();
gSystem->Sleep(200);
gSystem->ProcessEvents();
gSystem->DispatchOneEvent(true);
bool ok = askyn("Agree with the fit ? (Y/N)");
delete theApp;
return ok;
}
int main(int argnc, char **argv)
{
FILE *tabfile = nullptr;
cout<<" -- This program is used to create lookup tables for incub8r and 4play replay tasks"<<endl;
cout<<" Incub8r and 4play allow non-linear gains"<<endl;
cout<<" This is an adaptation of the lufwhm program of the Radware package"<<endl;
if(argnc != 2)
{
cout<<" -- To generate a binning file, use the command: ./gen_binning Conf_binning.txt"<<endl;
cout<<" --> see an example of conf file in conf/Conf_binning.txt"<<endl;
return 1;
}
vector< pair<float,float> > calibpoints;
pair<int, int> Range(0,8191);
Float_t NumberOfChannelPerFWHM = 1;
string line;
TString Buffer;
ifstream FileConf;
FileConf.open(argv[1], std::ifstream::in);
if(!FileConf) {
cout<<"No binning conf file ==> Exit"<<endl;
return 1;
}
int linenb=0;
while(FileConf) {
getline(FileConf,line);
linenb++;
Buffer = line;
TObjArray *loa=Buffer.ReplaceAll("\t"," ").Tokenize(" ");
if(Buffer.BeginsWith("#") || Buffer.Length()==0 ) continue;
if(Buffer.BeginsWith("TabFile")) {
if(loa->GetEntries() != 2) {
cout<<"Bad entry at line "<<linenb<<" ; skipped"<<endl; continue;
}
TabName = loa->At(1)->GetName();
if(!TabName.EndsWith(".tab")) TabName.Append(".tab");
}
if(Buffer.BeginsWith("FWHM")) {
if(loa->GetEntries() != 3) {
cout<<"Bad entry at line "<<linenb<<" ; skipped"<<endl; continue;
}
pair<float,float> value (((TString)loa->At(1)->GetName()).Atof(),((TString)loa->At(2)->GetName()).Atof());
calibpoints.push_back(value);
}
if(Buffer.BeginsWith("Precision")) {
if(loa->GetEntries() != 2) {
cout<<"Bad entry at line "<<linenb<<" ; skipped"<<endl; continue;
}
NumberOfChannelPerFWHM = ((TString)loa->At(1)->GetName()).Atof();
}
if(Buffer.BeginsWith("CompressionFactor")) {
if(loa->GetEntries() != 2) {
cout<<"Bad entry at line "<<linenb<<" ; skipped"<<endl; continue;
}
cfactor = ((TString)loa->At(1)->GetName()).Atof();
}
if(Buffer.BeginsWith("RangeMax")) {
if(loa->GetEntries() != 2) {
cout<<"Bad entry at line "<<linenb<<" ; skipped"<<endl; continue;
}
Range.first = 0;
Range.second = ((TString)loa->At(1)->GetName()).Atoi();
}
if(Buffer.BeginsWith("FixedBinning")) {
fFixedBinning = true;
}
delete loa;
}
tabfile = open_new_file(((char*)TabName.Data()), 1);
cout<<"Output tab file: "<<TabName.Data()<<endl;
chs = NumberOfChannelPerFWHM;
iclo = Range.first;
ichi = Range.second*cfactor;
cout<<" *-* "<<chs<<" channels per FWHM"<<endl;
cout<<" *-* Compression factor: "<<cfactor<<endl;
cout<<" *-* Channels Range: "<<iclo<<" -- "<<ichi<<endl;
if(fFixedBinning) {
cout<<" *-* Fixed binning used"<<endl;
wa2=1.;
wb2=0.;
wc2=0.;
}
else if(CalcFWHMParameters(calibpoints) == false)
return 1;
lo = iclo;
ich2 = 0;
while (lo <= ichi) {
++ich2;
x = ((float) lo + 0.5f) / 1e3f;
w = sqrt(wa2 + wb2 * x/cfactor + wc2 * x/cfactor * x/cfactor) * cfactor;
int incr = (int) ((w / chs) + 0.5);
if(incr==0) incr = 1;
lo += incr;
}
if (!askyn("No. of chs. in cube = %d --- Okay? (Y/N)", ich2))
return 1;
nclook = ichi+1;
for (j = 0; j < nclook; ++j) {
looktab[j] = 0;
}
lo = iclo;
for (i = 1; i < ich2+1; ++i) {
x = ((float) lo + 0.5f) / 1e3f;
w = sqrt(wa2 + wb2 * x/cfactor + wc2 * x/cfactor * x/cfactor) * cfactor;
hi = lo + (int) ((w / chs) + 0.5);
if (hi > nclook) hi = nclook;
for (j = lo; j < hi; ++j) {
looktab[j] = i;
}
wid_spec[i - 1] = (float) (hi - lo)/cfactor;
if ((lo = hi) == nclook) break;
}
/* write out look-up table to disk file */
lookmin = 1;
lookmax = ich2;
rewind(tabfile);
#define W(a,b) { memcpy(buf + rl, a, b); rl += b; }
W(&nclook, 4);
W(&lookmin, 4);
W(&lookmax, 4);
#undef W
if (put_file_rec(tabfile, buf, rl) ||
put_file_rec(tabfile, looktab, 2*nclook)) {
file_error("write to", ((char*)TabName.Data()));
return 1;
}
fclose(tabfile);
TString SpectrumFileName = TabName.Copy().Append(".spe");
tell("The \"Width Spectrum\" has %d channels, and contains the width (in ADC-chs) for each cube channel.\n", ich2);
wspec((char*)SpectrumFileName.Data(), wid_spec, ich2);
tell("=> Width Spectrum written as file %s\n", SpectrumFileName.Data());
return 0;
}
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
int
compress3d (unsigned int *in, unsigned char *out)
{
static int sum_limits[] = { 16,32,64,128,
256,512,1024,2048,
4096,8192,17000,35000,
80000,160000,320000,640000 };
static int bitmask[] = { 0,1,3,7,15,31,63,127,255,511,1023,2047,
4095,8191,16383,32767 };
int cube_length;
int sum=0, i, j;
int nplanes, noverflows, stmp, nbits;
unsigned char *bitptr, *ovrptr;