FSUDAQ_Qt6/analyzers/TEST.h
2024-06-04 11:59:18 -05:00

274 lines
9.1 KiB
C++

#ifndef TEST_h
#define TEST_h
/*********************************************
* This is online analyzer for RASIOR, ANL
*
* Created by Ryan @ 2023-10-16
*
* ******************************************/
#include "Analyser.h"
class TEST : public Analyzer{
public:
TEST(Digitizer ** digi, unsigned int nDigi, QMainWindow * parent = nullptr): Analyzer(digi, nDigi, parent){
SetUpdateTimeInSec(4.0);
RedefineEventBuilder({0}); // only builder for the 0-th digitizer.
tick2ns = digi[0]->GetTick2ns();
SetBackwardBuild(false, 100); // using normal building (acceding in time) or backward building, int the case of backward building, default events to be build is 100.
evtbder = GetEventBuilder();
evtbder->SetTimeWindow(500);
//========== use the influx from the Analyzer
influx = new InfluxDB("https://fsunuc.physics.fsu.edu/influx/");
dataBaseName = "testing";
SetUpCanvas(); // see below
};
void SetUpCanvas();
public slots:
void UpdateHistograms();
private:
MultiBuilder *evtbder;
Histogram2D * hPID;
Histogram2D * hXX; // X1 versus X2 : e[1] versus e[0]
Histogram2D * hYY; // Y1 versus Y2 : e[3] versus e[2]
Histogram1D * hXE; // X energy: e[0]+e[1]
Histogram1D * hYE; // Y energy: e[2]+e[3]
Histogram2D * hXYE; // 2D energy plot: e[2]+e[3] versus e[0]+e[1]
Histogram1D * hX; // X position:((e[0]-e[1])/(e[0]+e[1]))
Histogram1D * hY; // Y position:((e[2]-e[3])/(e[2]+e[3]))
Histogram2D * hXY; // 2D position plot: ((e[2]-e[3])/(e[2]+e[3])) versus ((e[0]-e[1])/(e[0]+e[1]))
Histogram2D * hXPE; // X position versus X energy: ((e[0]-e[1])/(e[0]+e[1])) versus (e[0]+e[1])
Histogram2D * hYPE; // Y position versus Y energy: ((e[2]-e[3])/(e[2]+e[3])) versus (e[2]+e[3])
//TH1F * hX1, * hX2, * hY1, *hY2;
Histogram2D * hXEdE1; //X energy versus dE signal 1
Histogram2D * hYEdE1; //Y energy versus dE signal 1
Histogram2D * hXEdE2; //X energy versus dE signal 2
Histogram2D * hYEdE2; //Y energy versus dE signal 2
Histogram1D * hX1, * hX2, * hY1, * hY2;
int chX1, chX2; // yellow, Red
int chY1, chY2; // Blue, White
int tick2ns;
float dE, E;
unsigned long long dE_t, E_t;
float e0, e1, e2, e3, dE1, dE2;
unsigned long long t0, t1, t2, t3, dE1_t, dE2_t;
//unsigned Int_t * energy;
//unsigned long energy;
/*
chX1 = 0; // left
chX2 = 1; // right
chY1 = 2; // top
chY2 = 3; // bottom
*/
};
inline void TEST::SetUpCanvas(){
setGeometry(0, 0, 1500, 2000);
//============ histograms
hPID = new Histogram2D("Test", "E", "dE", 100, 0, 11500, 100, 0, 11500, this);
layout->addWidget(hPID, 0, 0);
hXY = new Histogram2D("2D position plot", "X position", "Y position", 100, -100, 100, 100, -100, 100, this);
layout->addWidget(hXY, 0, 1);
hXX = new Histogram2D("X1 versus X2", "X2", "X1", 100, 0, 12500, 100, 0, 12500, this);
layout->addWidget(hXX, 0, 2);
hYY = new Histogram2D("Y1 versus Y2", "Y2", "Y1", 100, -2000, 12500, 100, -2000, 12500, this);
layout->addWidget(hYY, 0, 3);
hXE = new Histogram1D("X energy", "Ex", 300, -500, 25000, this);
layout->addWidget(hXE, 1, 0);
hYE = new Histogram1D("Y energy", "Ey", 300, -500, 25000, this);
layout->addWidget(hYE, 1, 1);
hX = new Histogram1D("X position", "X", 300, -1, 1, this);
layout->addWidget(hX, 1, 2);
hY = new Histogram1D("Y position", "Y", 300, -1, 1, this);
layout->addWidget(hY, 1, 3);
hXPE = new Histogram2D("X energy versus X position", "X position", "X energy", 100, -1, 1, 100, -2000, 25000, this);
layout->addWidget(hXPE, 0, 4);
hYPE = new Histogram2D("Y energy versus Y position", "Y position", "Y energy", 100, -1000, 1000, 100, -2000, 30000, this);
layout->addWidget(hYPE, 1, 4);
hXEdE1 = new Histogram2D("X energy versus dE signal 1", "Ex", "dE signal 1", 100, -500, 25000, 100, -200, 12500, this);
layout->addWidget(hXEdE1, 2, 0);
hYEdE1 = new Histogram2D("Y energy versus dE signal 1", "Ey", "dE signal 1", 100, -500, 25000, 100, -2000, 12500, this);
layout->addWidget(hYEdE1, 2, 1);
hXEdE2 = new Histogram2D("X energy versus dE signal 2", "Ex", "dE signal 2", 100, -500, 25000, 100, -500, 12500, this);
layout->addWidget(hXEdE2, 2, 2);
hYEdE2 = new Histogram2D("Y energy versus dE signal 2", "Ey", "dE signal 2", 100, -500, 25000, 100, -2000, 12500, this);
layout->addWidget(hYEdE2, 2, 3);
hXYE = new Histogram2D("2D energy plot", "Ex", "Ey", 100, 0, 25000, 100, 0, 30000, this);
layout->addWidget(hXYE, 2, 4);
}
inline void TEST::UpdateHistograms(){
if( this->isVisible() == false ) return;
BuildEvents(false); // call the event builder to build events
//============ Get events, and do analysis
long eventBuilt = evtbder->eventBuilt;
if( eventBuilt == 0 ) return;
//============ Get the cut list, if any
QList<QPolygonF> cutList = hPID->GetCutList();
const int nCut = cutList.count();
unsigned long long tMin[nCut] = {0xFFFFFFFFFFFFFFFF}, tMax[nCut] = {0};
unsigned int count[nCut]={0};
QList<QPolygonF> cutList1 = hXX->GetCutList();
const int nCut1 = cutList1.count();
unsigned long long tMin1[nCut1] = {0xFFFFFFFFFFFFFFFF}, tMax1[nCut1] = {0};
unsigned int count1[nCut1]={0};
//============ Processing data and fill histograms
long eventIndex = evtbder->eventIndex;
long eventStart = eventIndex - eventBuilt + 1;
if(eventStart < 0 ) eventStart += MaxNEvent;
for( long i = eventStart ; i <= eventIndex; i ++ ){
std::vector<Hit> event = evtbder->events[i];
//printf("-------------- %ld\n", i);
if( event.size() == 0 ) return;
for( int k = 0; k < (int) event.size(); k++ ){
//event[k].Print();
if( event[k].ch == 8 ) {dE = event[k].energy; dE_t = event[k].timestamp;} // Surface Barrier dE detector
if( event[k].ch == 8 ) {E = event[k].energy; E_t = event[k].timestamp;} // Surface Barrier E detector
if( event[k].ch == 8 ) {e0 = event[k].energy; t0 = event[k].timestamp;} //
if( event[k].ch == 8 ) {e1= event[k].energy; t1 = event[k].timestamp;} // The 4 output signals from the
if( event[k].ch == 8 ) {e2 = event[k].energy; t2 = event[k].timestamp;} // position sensitive E detector
if( event[k].ch == 8 ) {e3= event[k].energy; t3 = event[k].timestamp;} //
if( event[k].ch == 8 ) {dE1 = event[k].energy; dE1_t = event[k].timestamp;} // The 2 output signals from the
if( event[k].ch == 8 ) {dE2= event[k].energy; dE2_t = event[k].timestamp;} // square dE detector
}
// printf("(E, dE) = (%f, %f)\n", E, dE);
hPID->Fill(E + RandomGauss(0, 100), dE + RandomGauss(0, 100)); // x, y
//hXX->Fill(e1 + RandomGauss(0, 100), e0 + RandomGauss(0, 100)); //
hXX->Fill(e1, e0 ); //
hYY->Fill(e3 + RandomGauss(0, 100), e2 + RandomGauss(0, 100));
hXY->Fill(((e0-e1)/(e0+e1)) + RandomGauss(0, 100),((e2-e3)/(e2+e3)) + RandomGauss(0, 100));
hXE->Fill(e0+e1);
hYE->Fill(e2+e3);
hX->Fill(((e0-e1)/(e0+e1)));
hY->Fill(((e2-e3)/(e2+e3)));
hXPE->Fill(((e0-e1)/(e0+e1)) + RandomGauss(0, 100),(e0+e1) + RandomGauss(0, 100));
hYPE->Fill(((e2-e3)/(e2+e3)) + RandomGauss(0, 100),(e2+e3) + RandomGauss(0, 100));
hXEdE1->Fill((e0+e1)+ RandomGauss(0, 100),dE1 + RandomGauss(0, 100));
hYEdE1->Fill((e2+e3) + RandomGauss(0, 100),dE1 + RandomGauss(0, 100));
hXEdE2->Fill((e0+e1)+ RandomGauss(0, 100),dE2 + RandomGauss(0, 100));
hYEdE2->Fill((e2+e3)+ RandomGauss(0, 100),dE2 + + RandomGauss(0, 100));
hXYE->Fill((e0+e1) + RandomGauss(0, 100),(e2+e3) + RandomGauss(0, 100));
//check events inside any Graphical cut and extract the rate
for(int p = 0; p < cutList.count(); p++ ){
if( cutList[p].isEmpty() ) continue;
if( cutList[p].containsPoint(QPointF(E, dE), Qt::OddEvenFill) ){
if( dE_t < tMin[p] ) tMin[p] = dE_t;
if( dE_t > tMax[p] ) tMax[p] = dE_t;
count[p] ++;
//printf(".... %d \n", count[p]);
}
}
for(int p = 0; p < cutList1.count(); p++ ){
if( cutList1[p].isEmpty() ) continue;
if( cutList1[p].containsPoint(QPointF(e1, e0), Qt::OddEvenFill) ){
if( t0 < tMin1[p] ) tMin1[p] = t0;
if( t0 > tMax1[p] ) tMax1[p] = t0;
count1[p] ++;
}
}
}
for(int p = 0; p < cutList1.count(); p++ ){
printf("hXX.... %d %d \n", p, count1[p]);
}
hPID->UpdatePlot();
hXY->UpdatePlot();
hXX->UpdatePlot();
hYY->UpdatePlot();
hXE->UpdatePlot();
hYE->UpdatePlot();
hX->UpdatePlot();
hY->UpdatePlot();
hXPE->UpdatePlot();
hYPE->UpdatePlot();
hXEdE1->UpdatePlot();
hYEdE1->UpdatePlot();
hXEdE2->UpdatePlot();
hYEdE2->UpdatePlot();
hXYE->UpdatePlot();
//========== output to Influx
QList<QString> cutNameList = hPID->GetCutNameList();
for( int p = 0; p < cutList.count(); p ++){
if( cutList[p].isEmpty() ) continue;
double dT = (tMax[p]-tMin[p]) * tick2ns / 1e9; // tick to sec
double rate = count[p]*1.0/(dT);
//printf("%llu %llu, %f %d\n", tMin[p], tMax[p], dT, count[p]);
//printf("%10s | %d | %f Hz \n", cutNameList[p].toStdString().c_str(), count[p], rate);
influx->AddDataPoint("Cut,name=" + cutNameList[p].toStdString()+ " value=" + std::to_string(rate));
influx->WriteData(dataBaseName);
influx->ClearDataPointsBuffer();
}
}
#endif