Skip to content

Configuarion file and liear fit correction component. #10

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
arhayrap wants to merge 2 commits into
base: master
Choose a base branch
from
+108 −25
Open

Configuarion file and liear fit correction component. #10

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
ef39b99
Latest updates.
arhayrap Feb 2, 2024
7400087
Latest changes
arhayrap Apr 26, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 2 additions & 2 deletions config/LecroyScope_v12.config
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -9,5 +9,5 @@ z_DUT -15, -12, -9, -6, -3, 0, 3, 6, 9, 12
3 +. 0.05 0.15 0. 0. LP2+G40 0.
4 +. 0.05 0.15 0. 0. LP2+G40 0.
5 +. 0.05 0.15 0. 0. LP2+G40 0.
6 +. 0.05 0.15 0. 0. LP2+G40 0.
7 +. 0.05 0.15 0. 0. LP2+G40 0.
6 + 0.05 0.15 0. 0. Re+G40 0.
7 +. 0.05 0.15 0. 0. LP2+G40 0.
2 changes: 2 additions & 0 deletions include/DatAnalyzer.hh
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -9,6 +9,7 @@
// SYS includes
#include <sys/types.h>
#include <sys/stat.h>
#include <typeinfo>

// ROOT INCLUDES
#include "TFile.h"
Expand Down Expand Up @@ -70,6 +71,7 @@ class DatAnalyzer {
virtual void Analyze();

float GetPulseIntegral(float *a, float *t, unsigned int i_st, unsigned int i_stop); //returns charge in pC asssuming 50 Ohm termination
int GetMinIndex(float* array);
unsigned int GetIdxClosest(float value, float* v, unsigned int i_st, int direction=+1);
unsigned int GetIdxFirstCross(float value, float* v, unsigned int i_st, int direction=+1);
void AnalyticalPolinomialSolver(unsigned int Np, float* in_x, float* in_y, unsigned int deg, float* &out_coeff, float* err = 0);
Expand Down
15 changes: 15 additions & 0 deletions setup_cmslpc_python3.sh
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,15 @@
#!/bin/bash
echo Setting up CMS TimingDAQ.

source /cvmfs/sft.cern.ch/lcg/releases/LCG_103/Python/3.9.12/x86_64-ubuntu2004-gcc9-opt/Python-env.sh
source /cvmfs/sft.cern.ch/lcg/releases/LCG_103/ROOT/6.28.00/x86_64-ubuntu2004-gcc9-opt/ROOT-env.sh
source /cvmfs/sft.cern.ch/lcg/releases/LCG_103/scipy/1.8.0/x86_64-ubuntu2004-gcc9-opt/scipy-env.sh


# source /cvmfs/sft.cern.ch/lcg/releases/LCG_88/gcc/4.9.3/x86_64-slc6/setup.sh
# source /cvmfs/sft.cern.ch/lcg/releases/LCG_96python3/Python/3.6.5/x86_64-slc6-gcc8-opt/Python-env.sh
# source /cvmfs/sft.cern.ch/lcg/releases/LCG_88/ROOT/6.08.06/x86_64-slc6-gcc49-opt/bin/thisroot.sh

# source /cvmfs/sft.cern.ch/lcg/releases/LCG_96python3/Python/3.6.5/x86_64-slc6-gcc8-opt/Python-env.sh
# source /cvmfs/sft.cern.ch/lcg/releases/LCG_88/ROOT/6.08.06/x86_64-slc6-gcc49-opt/ROOT-env.sh
# source /cvmfs/sft.cern.ch/lcg/releases/LCG_88/scipy/0.18.1/x86_64-slc6-gcc49-opt/scipy-env.sh
112 changes: 89 additions & 23 deletions src/DatAnalyzer.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -263,11 +263,16 @@ void DatAnalyzer::Analyze(){
/*******************************
// -------------- Do linear fit
********************************/
if( config->channels[i].algorithm.Contains("Re") ) {
unsigned int i_min = GetIdxFirstCross(config->channels[i].re_bounds[0]*amp, channel[i], idx_min, -1);
unsigned int i_max = GetIdxFirstCross(config->channels[i].re_bounds[1]*amp, channel[i], i_min, +1);
if(config->channels[i].algorithm.Contains("Re") ) {
// cout<<sizeof(time[GetTimeIndex(i)])/sizeof(time[GetTimeIndex(i)][0])<<endl;
// unsigned int i_min = GetMinIndex(channel[i]); // GetIdxFirstCross(config->channels[i].re_bounds[0]*amp, channel[i], idx_min, -1); config->channels[i].re_bounds[0]*amp
unsigned int i_min = GetIdxFirstCross(config->channels[i].re_bounds[0]*amp, channel[i], idx_min, -1); // GetIdxClosest(min(channel[i]), channel[i], idx_min, -1);
unsigned int i_max = GetIdxFirstCross(config->channels[i].re_bounds[1]*amp, channel[i], i_min , +1);
float t_min = time[GetTimeIndex(i)][i_min];
float t_max = time[GetTimeIndex(i)][i_max];
// cout<<i_min<<" ----------- "<<i_max<<endl;
// cout<<t_min<<" ----------- "<<t_max<<endl;
// cout<<"==========================================="<<endl;

TF1* flinear = new TF1("flinear"+name, "[0]*x+[1]", t_min, t_max);
flinear->SetLineColor(2);
Expand All @@ -278,14 +283,61 @@ void DatAnalyzer::Analyze(){
pulse->Fit("flinear"+name, opt);
Re_slope = flinear->GetParameter(0);
Re_b = flinear->GetParameter(1);

// cout<<"==========================================="<<endl;
// cout<<config->channels[i].algorithm<<endl;
for ( auto f : config->constant_fraction ) {
var[Form("linear_RE_%d", (int)(100*f))][i] = (f*amp-Re_b)/Re_slope;
var[Form("linear_RE_%d", (int)(100*f))][i] = (f*amp-Re_b)/Re_slope + myTimeOffset;
}
// cout<<var[Form("linear_RE_%d", (int)(100*0.2))][i]<<endl;
// cout<<"==========================================="<<endl;

for ( auto thr : config->constant_threshold ) {
var[Form("linear_RE__%dmV", (int)(fabs(thr)))][i] = (thr-Re_b)/Re_slope;
var[Form("linear_RE__%dmV", (int)(fabs(thr)))][i] = (thr-Re_b)/Re_slope + myTimeOffset;
}
cout<<"==========================================="<<endl;

// if (i_evt > 0 && i_evt < 10) { // Plotting the clock and the fitted line.
TCanvas * c = new TCanvas(Form("fit_results%d", i_evt), Form("fit_results%d", i_evt));
// pulse->Draw(Form("channels[%d]:time", i), Form("event==%d", i_evt));
const int bins = NUM_SAMPLES; // round(sizeof(channel[i])/sizeof(int));
if (var[Form("linear_RE_%d", (int)(100*0.2))][i] == 0) {
cout<<"Clock time stamp is 0: "<<0.2*amp<<" "<<Re_b<<endl;
}
// cout<<t_min<<" "<<t_max<<" "<<bins<<endl;

float line[bins] = {};
float t = t_min;
float time_array[bins] = {};
// float clock_waveform[bins] = {};
for (int b = 0; b < bins; b++) {
t += b*(t_max - t_min)/bins;
time_array[b] = t;
line[b] = Re_slope * t + Re_b;
// clock_waveform[b] = channel[i][b]
// cout<<t<<" "<<line[b]<<endl;
}
// [GetTimeIndex(i)]
// for (int k =0; k<bins; k++) {
// cout<<k<<") "<<typeid(time[k]).name()<<" "<<typeid(channel[k]).name()<<endl;
// }
TGraph * clock_plot = new TGraph(bins, time[GetTimeIndex(i)], channel[i]);
TGraph * line_plot = new TGraph(bins, time_array, line);
cout<<"The size of the time: "<<sizeof(time[GetTimeIndex(i)])/sizeof(time[GetTimeIndex(i)][0])<<endl;
cout<<"The size of the channel: "<<sizeof(channel[i])/sizeof(channel[i][0])<<endl;
line_plot->SetLineColor(2);
line_plot->SetLineWidth(2);
clock_plot->SetLineColor(4);
clock_plot->SetLineWidth(2);

clock_plot->Draw("");
line_plot->Draw("same");

c->SaveAs(Form("/home/daq/ETROC2_Test_Stand/module_test_sw/analysis/FIT_PLOTS_TEST/fit_results%d.png", i_evt));
delete clock_plot;
delete line_plot;
delete c;

// }

delete flinear;
}
Expand Down Expand Up @@ -1288,15 +1340,14 @@ unsigned int DatAnalyzer::GetIdxClosest(float value, float* v, unsigned int i_st

unsigned int DatAnalyzer::GetIdxFirstCross(float value, float* v, unsigned int i_st, int direction) {
unsigned int idx_end = direction>0 ? NUM_SAMPLES-1 : 0;
bool rising = value > v[i_st]? true : false; // Check if the given max value is greater than the given waveform amplitude at a given start time (i_st).

bool rising = value > v[i_st]? true : false;

// if it is: the value of the waveform needs to rise and otherwise it doesn't need so the function will return the i_st.
unsigned int i = i_st;

while( i != idx_end ) {
if(rising && v[i] > value) break;
else if( !rising && v[i] < value) break;
i += direction;
while( (i != idx_end)) { // loop over the time bins
if(rising && v[i] > value) break; // check if the rising variable is true and if the waveform value is going higher than the given amplitude value. And stops the loop.
else if( !rising && v[i] < value) break; // check if the rising variable is false and if the waveform value is still lower than the given amplitude value. And stops the loop.
i += direction; // Otherwise it moves to the next time bin.
}

return i;
Expand Down Expand Up @@ -1437,6 +1488,21 @@ int DatAnalyzer::TimeOverThreshold(Interpolator *voltage, double tThresh, double

};

int DatAnalyzer::GetMinIndex(float* array) {
float min = 10000;
int index = 10000;
for (int i = 0; i < (sizeof(array) / sizeof(array[0])); i++) {
if (array[i] < min) {
min = array[i];
cout<<i<<", ";
index = i;
}
}
cout<<endl;
cout<<min<<endl;
return index;
};

float DatAnalyzer::FrequencySpectrum(double freq, double tMin, double tMax, int ich, int t_index){

const int range = 0; // extension of samples to be used beyond [tMin, tMax]
Expand All @@ -1460,21 +1526,21 @@ float DatAnalyzer::FrequencySpectrum(double freq, double tMin, double tMax, int

float DatAnalyzer::FrequencySpectrum(double freq, double tMin, double tMax, unsigned int n_samples, float* my_channel, float* my_time)
{
const int range = 0; // extension of samples to be used beyond [tMin, tMax]
double deltaT = (my_time[n_samples - 1] - my_time[0])/(double)(n_samples - 1); // sampling time interval
double fCut = 0.5/deltaT; // cut frequency = 0.5 * sampling frequency from WST
int n_min = floor(tMin/deltaT) - range; // first sample to use
int n_max = ceil(tMax/deltaT) + range; // last sample to use
n_min = std::max(n_min,0); // check low limit
n_max = std::min(n_max, (int)(n_samples - 1)); // check high limit
const int range = 0;
double deltaT = (my_time[n_samples - 1] - my_time[0])/(double)(n_samples - 1);
double fCut = 0.5/deltaT;
int n_min = floor(tMin/deltaT) - range;
int n_max = ceil(tMax/deltaT) + range;
n_min = std::max(n_min,0);
n_max = std::min(n_max, (int)(n_samples - 1));
int n_0 = (n_min + n_max)/2;

TComplex s(0.,0.); // Fourier transform at freq
TComplex I(0.,1.); // i
TComplex s(0.,0.);
TComplex I(0.,1.);

for(int n = n_min; n <= n_max; n++)
{
s += deltaT*(double)my_channel[n]*TComplex::Exp(-I*(2.*TMath::Pi()*freq*(n-n_0)*deltaT));//maybe don't need n_0 here, I think it will just add a phase to the fourier transform
s += deltaT*(double)my_channel[n]*TComplex::Exp(-I*(2.*TMath::Pi()*freq*(n-n_0)*deltaT));
}
return s.Rho();
};