A double-sided silicon strip detector (DSSD) is used in nuclear physics to pinpoint the position and energy of incident particles. With orthogonal strip arrays on its front (x-strips) and back (y-strips), it detects charge carriers generated by a particle’s energy deposition. These carriers induce equal but opposite charges on the x- and y-strips, and amplifiers of opposite polarity produce positive signals of equal amplitude, reflecting the deposited energy.
For a particle hitting the intersection of front strip $x_i$ and back strip $y_j$, depositing energy $E$, the calibrated energy signals from both sides should match:
where $E_{xi}$ and $E_{yj}$ are the energies measured by strips $x_i$ and $y_j$.
The raw signal amplitudes (channel numbers) from strips $x_i$ and $y_j$, denoted $a_{xi}$ and $a_{yj}$, relate to the deposited energy via linear calibration:
$$ E_{xi} = g_{xi} a_{xi} + o_{xi}, \tag{1} $$$$ E_{yj} = g_{yj} a_{yj} + o_{yj}, \tag{2} $$where $g_{xi}, g_{yj}$ are gain factors and $o_{xi}, o_{yj}$ are offsets. Since the deposited energy is the same ($E = E_{xi} = E_{yj}$), equating (1) and (2) gives a linear relationship between amplitudes:
$$ a_{xi} = k_{yj} a_{yj} + b_{yj}, \tag{3} $$where: $$ k_{yj} = \frac{g_{yj}}{g_{xi}}, \quad b_{yj} = \frac{o_{yj} - o_{xi}}{g_{xi}}. $$
Here, $k_{yj}$ and $b_{yj}$ are normalization coefficients that scale the amplitude $a_{yj}$ (from the "calibrated" strip $y_j$) to match $a_{xi}$ (from the "reference" strip $x_i$). A 2D plot of $a_{xi}$ vs. $a_{yj}$ shows a linear trend, and a linear fit yields $k_{yj}$ and $b_{yj}$. This calibration ensures consistent energy measurements across strips, crucial for accurate position and energy reconstruction.
When two particles strike the DSSD simultaneously at positions $(x_1, y_1)$ and $(x_2, y_2)$ with distinct energies $e_1 = e_{x1} = e_{y1}$ and $e_2 = e_{x2} = e_{y2}$ ($e_1 \neq e_2$), signals are recorded from x-strips $x_1, x_2$ and y-strips $y_1, y_2$. This leads to ambiguity, as possible pairings are:
The correct pairing produces the smallest energy differences. After normalizing amplitudes using equation (3), the proper combination of $(x_1, y_1)$ and $(x_2, y_2)$ minimizes $|e_x - e_y|$ among all x- and y-strip pairings, enabling precise position reconstruction.
Reaction: 25 MeV/A $^{16}$C + $^{9}$Be, measuring the multiple charged particle products from the fragmentation of $^{16}$C on the target.
Detector Setup: A telescope array of 3 DSSDs + CsI, positioned at zero angle behind the target, facing the beam direction.
D1: 32 x 32 - D2: 32 x 32 - D3: 32 x 32
SSD: 3 pieces
CsI: 2 x 2
Geometry:
Target: Diameter: ϕ=30\phi = 30\phi = 30mm; z = 0 mm
D1: z = 156 mm
D2: z = 176 mm
D3: z = 196 mm
Trigger: Requires the x-side of both D1 and D2 to have signals from two or more strips
Event Preprocessing: Pedestal events have been removed.
// xenergy, yenergy, xtime
d1x[32], d1y[32], d1t[32];
d2x[32], d2y[32], d2t[32];
d3x[32], d3y[32], d3t[32];
//hit, energy, time, strip
d1xhit, d1xe[d1xhit], d1xt[d1xhit], d1xs[d1xhit];
d2xhit, d2xe[d2xhit], d2xt[d2xhit], d2xs[d2xhit]
d3xhit, d3xe[d3xhit], d3xt[d3xhit], d3xs[d3xhit]
d1yhit, d1ye[d1yhit], d1yt[d1yhit], d1ys[d1yhit];
d2yhit, d2ye[d2yhit], d2yt[d2yhit], d2ys[d2yhit]
d3yhit, d3ye[d3yhit], d3yt[d3yhit], d3ys[d3yhit]
//%jsroot on
TFile *ipf = new TFile("./data/data_16C.root");
TTree *tree = (TTree*)ipf->Get("tree");
TCanvas *c1=new TCanvas("c1","c1");
Sorted in descending order of energy.
Mostly adjacent strip events.
tree->Scan("d1xe:d1xs:d1ye:d1ys","","",10,1);
*********************************************************************** * Row * Instance * d1xe * d1xs * d1ye * d1ys * *********************************************************************** * 1 * 0 * 4800 * 17 * 5775 * 18 * * 1 * 1 * 982 * 18 * * * * 2 * 0 * 3069 * 23 * 3352 * 11 * * 2 * 1 * 341 * 24 * 111 * 10 * * 3 * 0 * 5325 * 23 * 4522 * 15 * * 3 * 1 * * * 922 * 16 * * 4 * 0 * 5822 * 21 * 5688 * 15 * * 4 * 1 * * * 271 * 14 * * 5 * 0 * 4048 * 11 * 5213 * 12 * * 5 * 1 * 1149 * 12 * * * * 6 * 0 * 3719 * 14 * 4255 * 16 * * 6 * 1 * 374 * 13 * * * * 7 * 0 * 4546 * 12 * 4961 * 9 * * 7 * 1 * 297 * 11 * * * * 8 * 0 * 4051 * 22 * 3877 * 9 * * 8 * 1 * * * 329 * 10 * * 9 * 0 * 6720 * 13 * 6513 * 12 * * 9 * 1 * 142 * 10 * 418 * 13 * * 9 * 2 * 109 * 14 * * * * 10 * 0 * 3873 * 20 * 3007 * 8 * * 10 * 1 * * * 974 * 9 * ***********************************************************************
The majority of events on both the x and y sides are characterized by twofold multiplicity.
A small number of events exhibit less than twofold multiplicity on both sides, as events with energy below the pedestal threshold have been filtered out.
tree->Draw("d1xhit:d1yhit>>(15,0,15,15,0,15)","","colz");
gPad->SetLogz();
c1->Draw();
A clear linear correlation region is observed in the central area.
Scattered points outside the linear region result from energy sharing between adjacent strips on the x/y sides.
gPad->SetLogy(0);
tree->Draw("d1x[12]:d1y[13]>>(1000,0,8000,1000,0,8000)","","colz");
c1->Draw();//front-back correlation
To mitigate the impact of energy sharing between adjacent strips, two approaches are used:
Single Multiplicity Requirement: Require exactly one hit (single multiplicity) on both x and y sides. This significantly reduces scattered points but results in lower event statistics due to trigger conditions.
No Multiplicity Restriction: Instead of limiting multiplicity, require that the selected strips on the x and y sides have no energy sharing with adjacent strips.
For d1x[12]: Ensure strips 11 and 13 have no signal sharing.
For d1y[13]: Ensure strips 12 and 14 have no signal sharing.
Applying these improved correlation conditions significantly increases the event statistics.
c1->Clear();
TCut chit1="d1xs==12 && d1ys==13 && d1xhit==1 && d1yhit==1";
tree->Draw("d1xe:d1ye>>(1000,0,8000,1000,0,8000)",chit1,"");
c1->Draw();
c1->Clear();
TCut cveto = "d1x[11]<50 && d1x[13]<50 && d1y[12]<50 && d1y[14]<50";//no crosstalk from neighbouring strips
TCut c1213 = "d1x[12]>200 &&d1y[13]>200&& d1y[13]<8000" && cveto;
tree->Draw("d1x[12]:d1y[13]>>h2(1000,0,8000,1000,0,8000)", c1213, "");//front-back correlation with crosstalk rejection
c1->Draw();
TGraph
with events from a TTree
using the sequence of:TTree::Draw("var1:var2", "cuts", "goff")
TGraph *gr = new TGraph(tree->GetSelectedRows(), tree->GetV2(), tree->GetV1())
.TGraph *gr = new TGraph(tree->GetSelectedRows(), tree->GetV2(), tree->GetV1());
gr->SetMarkerSize(0.2);
gr->Draw("A*");//draw point
c1->Draw();
Fit Results Deviate from the Actual X-Y Correlation Line
A significant number of outliers exist outside the linear correlation, which are not caused by statistical fluctuations but by other mechanisms, adversely affecting the fit results.
The fit can be improved by manually selecting a fitting region using TCutG
.
TF1 *fp1 = new TF1("fp1", "pol1", 200, 8000);
gr->Fit(fp1,"R");
gr->Draw("A*");
c1->Draw();
**************************************** Minimizer is Linear / Migrad Chi2 = 8.70489e+08 NDf = 586 p0 = 638.618 +/- 70.9968 p1 = 0.516471 +/- 0.0314706
In the presence of outliers that do not come from the same data-generating process as the rest of the data, least squares estimation is inefficient and can be biased. Robust regression methods are designed to be not overly affected by violations of assumptions by the underlying data-generating process. common situation in which robust estimation is used occurs when the data contain outliers. https://root.cern/doc/master/fitLinearRobust_8C.html
ROOT only supports robust fitting for polynomial functions.
Double_t p0,p1,p2;
c1->Clear();
gr->Fit(fp1, "ROB");
//The option "rob=0.75" means that we want to use robust fitting and
//we know that at least 75% of data is good points (at least 50% of points
//should be good to use this algorithm). If you don't specify any number
//and just use "rob" for the option, default value of (npoints+nparameters+1)/2
//will be taken
p0=fp1->GetParameter(0);
p1=fp1->GetParameter(1);
gr->Draw("A*");
c1->Draw();
**************************************** Minimizer is Linear / Robust Chi2 = 1.20668e+09 NDf = 586 p0 = 12.9492 p1 = 0.983987
The linear function's fit is less satisfactory in areas with larger x values.
TString stree;
stree.Form("d1y[13]:d1x[12]-(%lf*d1y[13]+%lf)>>ha(100,-40,40,1000,0,8000)", p1, p0);//%lf-Double_t
tree->Draw(stree.Data(), c1213, "");
c1->Draw();
TF1 *fp2 = new TF1("fp2","pol2",200,8000);
gr->Fit(fp2,"R ROB");
p0 = fp2->GetParameter(0);
p1 = fp2->GetParameter(1);
p2 = fp2->GetParameter(2);
stree.Form("d1y[13]:d1x[12]-(%e*d1y[13]*d1y[13]+%lf*d1y[13]+%f)>>ha(100,-40,40,1000,0,8000)", p2, p1, p0);
tree->Draw(stree.Data(), c1213, "");
c1->Draw();
**************************************** Minimizer is Linear / Robust Chi2 = 1.20693e+09 NDf = 585 p0 = 10.6689 p1 = 0.987026 p2 = -5.43837e-07
Save y~x events from DSSD1/2/3, excluding the influence of adjacent strips, into d1/2/3xy.root files:
//load a cut file.
gROOT->Macro("./data/cut1.C");
TCutG *cut1 = (TCutG *)gROOT->GetListOfSpecials()->FindObject("cut1");
//Checks if a pixel meets isolation criteria (no adjacent signal sharing)
bool IsIsolatedPixel(int dxe[32], int dye[32], int ix, int iy) {
...
}
//input ROOT file
TFile *fin = new TFile("data_16C.root");
TTree *tree =(TTree *)fin->Get("tree");
Int_t d1x[32], d2x[32], d3x[32];
Int_t d1y[32], d2y[32], d3y[32];
...
//output ROOT file
int ix, iy;
Int_t xe, ye;
TFile *fout = new TFile("d1dxy.root","recreate");
...
for (Long64_t jentry = 0; jentry < nentries; jentry++) {
tree->GetEntry(jentry);
bool bcut;
for (int i = 0; i < 32; i++) {
for (int j = 0; j < 32; j++) {
bcut = cut1->IsInside(d1x[i], d1y[j]);
bcut *= IsIsolatedPixel(d1x,d1y,i,j);
if (!bcut) continue;
...
tout->Fill();
}
}
if (jentry % 10000 == 0 ) cout << ".";
}
...
TCanvas *c1 = new TCanvas("c1","c1");
TFile *fin = new TFile("./data/d1xy.root");
TTree *tree =(TTree *)fin->Get("tree");
Warning in <TCanvas::Constructor>: Deleting canvas with same name: c1
tree->Draw("ye:xe>>(4000,0,8000,4000,0,8000)","","colz");
c1->SetLogz();
c1->Draw();
tree->Scan("iy:ix:ye:xe","","",10,1);
************************************************************ * Row * iy * ix * ye * xe * ************************************************************ * 1 * 23 * 4 * 501 * 629 * * 2 * 7 * 10 * 611 * 471 * * 3 * 23 * 10 * 501 * 471 * * 4 * 18 * 3 * 1171 * 1178 * * 5 * 14 * 8 * 507 * 498 * * 6 * 23 * 26 * 1387 * 1351 * * 7 * 21 * 29 * 459 * 444 * * 8 * 20 * 25 * 445 * 491 * * 9 * 27 * 25 * 484 * 491 * * 10 * 20 * 27 * 445 * 437 * ************************************************************