ROOT Tutorial II：TTree、关联分析与事件选择
建议学完第一章后阅读。编程基础见 Tutorial I。
以 ΔE–E 望远镜为例，将同一粒子的探测器信号保存在 TTree 中，再通过逻辑 cut 和 graphical cut 选择粒子、构建能谱与 PID 投影。
按顺序运行 notebook。两种语言使用相同数据模型和筛选条件，选择一种即可。
0. 在 Jupyter 中启动 ROOT C++¶
选择 ROOT C++ kernel，按顺序从头运行 notebook。后面的部分单元格会读取前面生成的文件；重启 kernel 会清除内存中的 C++ 对象，但不会删除磁盘上的 ROOT 文件。
#include <iostream>
#include <cmath>
#include <algorithm>
#include <TROOT.h>
#include <TCanvas.h>
#include <TFile.h>
#include <TTree.h>
#include <TRandom3.h>
#include <TH1.h>
#include <TH2.h>
#include <TCut.h>
#include <TCutG.h>
#include <TLegend.h>
#include <TLatex.h>
#include <TStyle.h>
std::cout << "ROOT version: " << gROOT->GetVersion() << std::endl;
ROOT version: 6.40.02
头文件声明本页用到的 ROOT 类。TCut 保存逻辑条件，TCutG 保存二维多边形区域；后面分别演示。
%jsroot on 是 Jupyter 命令，用于开启交互画布，不是 C++ 语句。在 ROOT macro 中运行时去掉这一行。
%jsroot on
gStyle->SetOptStat(0);
TCanvas *c1 = new TCanvas("c1", "Telescope analysis", 800, 520);
各图沿用这张画布。SetOptStat(0) 隐藏默认统计框，计数在需要时直接输出。
1. 探测器背景与事例记录
ΔE–E 望远镜测量什么
带电粒子进入硅后，通过电离和激发损失动能。硅探测器收集电离产生的电荷，信号经能量刻度后给出粒子在该层的能量沉积。将两层探测器前后排列，前层测量能量损失 ΔE，后层测量剩余能量 E，这就是 ΔE–E 望远镜。
本例中的粒子穿透前层、停止在后层：穿透表示离开该层时仍有动能，停止表示在该层耗尽剩余动能。忽略层间损失时，入射动能满足 \(E_0=\Delta E+E\)。如果粒子也穿出后层，两层信号之和就不再给出完整的入射能量。
为什么不同同位素形成条带
在非相对论能区，忽略缓慢变化的项时，阻止本领有 \(-dE/dx\propto Z^2/v^2\) 的变化趋势。这里 Z 和 A 表示入射粒子的电荷数和质量数。由 \(v^2\propto E_0/A\)，薄层中的能损近似随 \(Z^2A/E_0\) 变化。因此同一种粒子的 ΔE 和 E 彼此关联，不同种类形成不同的弯曲条带。相同 Z、相同速度的同位素阻止本领近似相同；在相同总动能下，它们的速度不同，能损也不同。
只有 ΔE 的一维分布会把不同能量、不同粒子的信号叠在一起。保留同一粒子的 E 后，才能在二维图中沿条带选择某种粒子。条带有宽度，来自能损涨落、测量分辨等因素，邻近条带也可能重叠。
event、entry 与 branch
一次入射粒子的完整记录称为事例（event）。TTree 通常用一个 entry 保存一个事例，用不同 branch 保存字段。两层信号放在同一个数组 branch 的两个元素中；读取某个 entry 时，两者仍属于同一粒子。直方图只累积 bin 计数，TTree 则保留逐事例数据，方便以后改变 cut 或重新计算观测量。
本例中的望远镜事例¶
设 p、d、t、³He、⁴He、⁶He 混合轻离子入射到薄硅 ΔE 探测器，后接足够厚、能够阻止粒子的 E 探测器。第一层测得 signal[0]，第二层测得剩余能量 signal[1]。两个值来自同一个粒子，因此保存在同一个 entry 中。
每条模拟记录包含：
eventID：按顺序生成的事例编号；A、Z：生成粒子的质量数和电荷数；E0：生成的入射动能，单位为 MeV；signal[0]：薄探测器测得的能量损失，单位为 MeV；signal[1]：厚探测器测得的能量，单位为 MeV。
生成的 (A, Z) = (1, 1)、(2, 1)、(3, 1)、(3, 2)、(4, 2)、(6, 2) 依次对应 p、d、t、³He、⁴He、⁶He。A、Z 和 E0 是模拟真值（simulation truth），可用于核对示例；普通实验数据不会自动给出这些真值。后面先根据测量信号筛选，再用真值标签检查。
生成示意性的望远镜数据
为练习 TTree，先生成一组简单的双层探测器数据。根据非相对论近似下 \(dE/dx\propto Z^2/v^2\) 和 \(v^2\propto E_0/A\) 的变化趋势，采用示意响应
$$\Delta E_{\rm true}=K\frac{Z^2 A}{E_0},\qquad E_{\rm true}=E_0-\Delta E_{\rm true},$$
其中 \(K=8\) MeV²，\(20\le E_0\le80\) MeV。假设粒子停止在第二层，再分别给两个信号加入 Gaussian 测量涨落。下面的参数仅供示例使用，不是实际探测器的刻度结果。这个解析模型只用于生成便于练习的条带，不代替第一章基于射程数据的能损计算；作业 4.1 使用作业 1.2 的三层探测器能量沉积。
2. 创建并填充 TTree
先打开输出文件，定义 branch，再逐事例填充。每调用一次 Fill()，就把当前全部字段保存为一条记录。
TFile *fout = new TFile("tree_telescope_cpp.root", "RECREATE");
TTree *tree_out = new TTree("telescope", "two-detector charged-particle telescope events");
Int_t eventID;
Int_t A;
Int_t Z;
Float_t E0;
Float_t signal_buffer[2];
TFile(filename, mode) 创建或打开 ROOT 文件。RECREATE 创建新文件，并覆盖已有同名文件，因此仅用于允许覆盖的输出文件。
TTree(name, title) 创建空 tree。telescope 是之后读取时使用的内部名称，title 用于说明内容。
这里的 C++ 变量作为内存 buffer。Int_t 是 ROOT 的 32 位有符号整数类型，Float_t 是 32 位浮点类型；signal_buffer[2] 为两层信号预留两个连续的 Float_t 元素。声明变量本身不会创建 branch 或保存 entry。
tree_out->Branch("eventID", &eventID, "eventID/I");
tree_out->Branch("A", &A, "A/I");
tree_out->Branch("Z", &Z, "Z/I");
tree_out->Branch("E0", &E0, "E0/F");
tree_out->Branch("signal", signal_buffer, "signal[2]/F");
Branch(name, address, leaflist) 的三个参数依次为 branch 名称、内存地址和数据结构。"A/I" 表示一个 Int_t；"signal[2]/F" 表示两个 Float_t。标量用 &A 传地址，数组名 signal_buffer 已表示首元素地址。内存类型和数组长度应与 leaf list 一致。
TRandom3 rng(2026);
const Int_t n_events = 60000;
const Double_t K_loss = 8.0; // 示意模型中的常数，单位为 MeV^2
const Int_t species_A[6] = {1, 2, 3, 3, 4, 6};
const Int_t species_Z[6] = {1, 1, 1, 2, 2, 2};
对每个粒子，先更新粒子标签、入射能量和两个信号，再调用一次 Fill()，将这些值共同保存为一条记录。
for (eventID = 0; eventID < n_events; ++eventID) {
Int_t species_index = rng.Integer(6);
A = species_A[species_index];
Z = species_Z[species_index];
E0 = rng.Uniform(20.0, 80.0);
Double_t dE_true = K_loss * Z * Z * A / E0;
Double_t E_true = E0 - dE_true;
Double_t sigma_dE = 0.03 + 0.04 * std::sqrt(dE_true);
Double_t sigma_E = 0.08 + 0.015 * std::sqrt(E_true);
signal_buffer[0] = std::max(0.0, rng.Gaus(dE_true, sigma_dE));
signal_buffer[1] = std::max(0.0, rng.Gaus(E_true, sigma_E));
tree_out->Fill();
}
std::cout << "Entries held in memory: "
<< tree_out->GetEntries() << std::endl;
Entries held in memory: 60000
Integer(6) 等概率选取六种同位素之一；Gaus(mean, sigma) 生成测量信号，其中两种 sigma 表达式都是本示意模型的假设。Fill() 将当前 branch buffer 复制到下一条记录；不要重新创建 buffer，只更新其中的值。
fout->cd();
tree_out->Write();
fout->Close();
std::cout << "Created tree_telescope_cpp.root" << std::endl;
Created tree_telescope_cpp.root
fout->cd() 将输出文件设为当前 ROOT 目录，Write() 保存 TTree 的 branch 结构和已填入的记录，Close() 完成文件写入并关闭。仅调用 Fill() 不能保证 TTree 已完整保存到磁盘，还需要完成写入步骤。
TTree 在 fout 为当前目录时创建，通常由该文件管理。关闭 fout 后，不要继续使用 tree_out；应像下面那样重新打开文件并取得新的 tree 指针。
3. 重新打开并查看 TTree
先检查文件对象、branch 名称与类型、数组长度和 entry 数，再进行分析。
TFile *fin = TFile::Open("tree_telescope_cpp.root", "READ");
fin->ls();
TTree *tree_in = nullptr;
fin->GetObject("telescope", tree_in);
TFile** tree_telescope_cpp.root TFile* tree_telescope_cpp.root KEY: TTree telescope;1 two-detector charged-particle telescope events
ls() 列出文件中的对象。按保存的名称 telescope 取出 TTree，并在分析期间保持输入文件打开。若文件或 tree 不存在，先核对文件名及对象列表。
tree_in->Print();
std::cout << "Number of entries: "
<< tree_in->GetEntries() << std::endl;
****************************************************************************** *Tree :telescope : two-detector charged-particle telescope events * *Entries : 60000 : Total = 1446637 bytes File Size = 806256 * * : : Tree compression factor = 1.79 * ****************************************************************************** *Br 0 :eventID : eventID/I * *Entries : 60000 : Total Size= 241125 bytes File Size = 84525 * *Baskets : 8 : Basket Size= 32000 bytes Compression= 2.85 * *............................................................................* *Br 1 :A : A/I * *Entries : 60000 : Total Size= 241053 bytes File Size = 39081 * *Baskets : 8 : Basket Size= 32000 bytes Compression= 6.16 * *............................................................................* *Br 2 :Z : Z/I * *Entries : 60000 : Total Size= 241053 bytes File Size = 25163 * *Baskets : 8 : Basket Size= 32000 bytes Compression= 9.56 * *............................................................................* *Br 3 :E0 : E0/F * *Entries : 60000 : Total Size= 241065 bytes File Size = 212170 * *Baskets : 8 : Basket Size= 32000 bytes Compression= 1.13 * *............................................................................* *Br 4 :signal : signal[2]/F * *Entries : 60000 : Total Size= 481883 bytes File Size = 444124 * *Baskets : 16 : Basket Size= 32000 bytes Compression= 1.08 * *............................................................................* Number of entries: 60000
Print() 显示 entry 数、branch 名称、类型和维度。这里 signal[2]/F 表示每条记录包含两个 Float_t 值，GetEntries() 返回 entry 总数。
tree_in->Show(10);
======> EVENT:10
eventID = 10
A = 2
Z = 1
E0 = 71.276
signal = 0.266728,
71.1437
Show(entry_number) 输出某条记录的所有 branch。entry 编号从 0 开始，因此 Show(10) 显示的是第十一条记录。它用于检查数据，不是提取分析结果的方法。
tree_in->Scan(
"eventID:A:Z:E0:signal[0]:signal[1]",
"",
"",
10,
0
);
************************************************************************************ * Row * eventID * A * Z * E0 * signal[0] * signal[1] * ************************************************************************************ * 0 * 0 * 2 * 1 * 65.996543 * 0.2959895 * 65.936271 * * 1 * 1 * 6 * 2 * 25.333940 * 7.5492363 * 17.606323 * * 2 * 2 * 6 * 2 * 57.311428 * 3.4103038 * 54.082107 * * 3 * 3 * 2 * 1 * 53.495006 * 0.3692665 * 53.215637 * * 4 * 4 * 4 * 2 * 38.148593 * 3.3014807 * 35.070316 * * 5 * 5 * 2 * 1 * 52.400577 * 0.3346374 * 52.181098 * * 6 * 6 * 1 * 1 * 69.937141 * 0.1737575 * 69.656311 * * 7 * 7 * 2 * 1 * 41.822879 * 0.4038591 * 41.643936 * * 8 * 8 * 1 * 1 * 21.835584 * 0.2760921 * 21.359891 * * 9 * 9 * 2 * 1 * 24.489458 * 0.6236841 * 23.949422 * ************************************************************************************
这里使用 Scan(expression, selection, option, number, first_entry)。用冒号分隔要显示的表达式；selection 和 option 留空表示不筛选、采用默认格式。最后两个参数指定从编号 0 开始显示十条记录。Scan() 适合抽查少量数据，不宜靠打印成千上万条记录完成分析或验证。
4. 二维关联、逻辑 cut 与 graphical cut
4.1 从二维关联图开始
Draw(expression, selection, option) 遍历 TTree。二维表达式写作 y:x，所以 signal[0]:signal[1] 表示横轴 E、纵轴 ΔE；TH2F 构造函数则先给 x 轴、再给 y 轴的 bin 设置。
TH2F *hPID = new TH2F("hPID", "Telescope;E (MeV);#DeltaE (MeV)",
850, 0, 85, 1100, 0, 11);
tree_in->Draw("signal[0]:signal[1]>>hPID", "", "COLZ");
c1->Draw();
>>hPID 将计算值填入指定直方图，COLZ 用颜色表示计数。这里 bin 宽度为 0.1 MeV × 0.01 MeV。直接重复 Draw 会先清空同名直方图；>>+hPID 才是累加。
氢同位素的 ΔE 较小，在全范围图中挤在底部。下面只放大坐标显示，不重新生成或筛选数据。
hPID->GetXaxis()->SetRangeUser(18, 38);
hPID->GetYaxis()->SetRangeUser(0, 1.5);
hPID->Draw("COLZ");
TLatex *labels = new TLatex();
labels->SetTextSize(0.035);
labels->DrawLatex(22, 0.22, "p");
labels->DrawLatex(22, 0.57, "d");
labels->DrawLatex(22, 0.92, "t");
c1->Draw();
局部图中由下到上是 p、d、t。增加 bin、放大坐标有助于看清原有结构，但不会改善探测器分辨；能量较高处的氢条带仍可能重叠。先用 SetRange(0, 0) 恢复全轴范围，后面继续分析全部事例。
hPID->GetXaxis()->SetRange(0, 0);
hPID->GetYaxis()->SetRange(0, 0);
tree_in->SetAlias("Etot", "signal[0]+signal[1]");
SetAlias 为表达式取别名。Etot 可用于后面的 Draw 和 cut，但不会增加一个 branch。
4.2 用 AND、OR、NOT 组合条件
cut 判断每条记录是否被接受。下面选择总能量在 20–40 MeV 或 60–80 MeV，并且两层信号均为正的事例。
| 写法 | 含义 | 例子 |
|---|---|---|
&&（AND） | 两项都成立 | Etot>=20 && Etot<40 |
||（OR） | 至少一项成立 | (Etot<40) || (Etot>=60) |
!（NOT） | 条件取反 | !(Etot>=40 && Etot<60) |
==、!= | 相等、不相等 | A==3 && Z==2（模拟中的 ³He 标签） |
无论外层代码是 PyROOT 还是 C++，传给 ROOT 的字符串表达式都采用上表写法。区间要写成两个比较，不写 20<Etot<40；混合 AND 和 OR 时，用括号明确组合。
TCut valid = "signal[0]>0 && signal[1]>0";
TCut low = "Etot>=20 && Etot<40";
TCut high = "Etot>=60 && Etot<80";
TCut windows = valid && (low || high);
std::cout << "Low window: " << tree_in->GetEntries(valid && low) << std::endl;
std::cout << "High window: " << tree_in->GetEntries(valid && high) << std::endl;
std::cout << "Either window: " << tree_in->GetEntries(windows) << std::endl;
std::cout << "Outside both: " << tree_in->GetEntries(valid && !(low || high)) << std::endl;
Low window: 20133 High window: 19885 Either window: 40018 Outside both: 19959
PyROOT 中可用字符串变量组合表达式；C++ 中的 TCut 是对条件字符串的封装，可用 &&、|| 组合。GetEntries(cut) 只统计满足条件的 entry，不画图。这两个能区不重叠，因此 OR 的计数等于两者相加；一般的 OR 不会重复计数同时满足两项的事例。
TH2F *hWindows = new TH2F("hWindows", "Two total-energy windows;E (MeV);#DeltaE (MeV)",
850, 0, 85, 1100, 0, 11);
tree_in->Draw("signal[0]:signal[1]>>hWindows", windows, "COLZ");
c1->Draw();
图中保留两段总能量区间，边界由 ΔE + E 决定，因此不是竖直线。筛选字段可以不同于绘图字段：例如三层探测器中，可用第三层 e[2] 的信号筛选，再画前两层 e[0]:e[1] 的关联。
4.3 用 TCutG 圈选 ³He 条带
能量区间适合用逻辑 cut 表达；沿弯曲条带选择粒子时，封闭多边形更直接，称为 graphical cut。下面沿 ³He 条带画一个有限范围的区域，再将区域内的事例用于其他分析。
TCutG(name, n) 创建 n 个顶点的 cut，SetPoint(i, x, y) 按轮廓顺序设置顶点，最后一个点与第一个点相同以闭合。SetVarX、SetVarY 说明横纵坐标对应哪个 tree 表达式；这里分别为后层 E 和前层 ΔE，不要与 Draw 的 y:x 顺序混淆。
Double_t xgate[] = {22, 30, 40, 50, 60, 70, 70, 60, 50, 40, 30, 22, 22};
Double_t ygate[] = {4.03, 3.24, 2.56, 2.12, 1.80, 1.55, 1.12, 1.31, 1.57, 1.98, 2.56, 3.43, 4.03};
TCutG *he3 = new TCutG("he3_cut", 13);
he3->SetVarX("signal[1]");
he3->SetVarY("signal[0]");
for (Int_t i = 0; i < 13; ++i) {
he3->SetPoint(i, xgate[i], ygate[i]);
}
hPID->Draw("COLZ");
he3->SetLineColor(kRed + 1);
he3->SetLineWidth(2);
he3->Draw("L SAME");
labels->DrawLatex(34, 3.7, "^{3}He gate");
c1->Draw();
顶点坐标来自这张图中的条带，只适用于本例。在 ROOT 桌面画布中，也可用图形编辑器的 CutG 工具逐点画出轮廓，再将默认名称 CUTG 改成有意义的名称。不同 notebook 前端的编辑功能不完全相同，直接设置顶点的写法便于复现和修改。
将 cut 名称 "he3_cut" 放入 selection，即可选择区域内的事例；也能继续写 "he3_cut && Etot<50" 或 "!he3_cut"。仅把轮廓画在图上不会筛选数据。
TH2F *hHe3 = new TH2F("hHe3", "Selected by the graphical cut;E (MeV);#DeltaE (MeV)",
850, 0, 85, 1100, 0, 11);
tree_in->Draw("signal[0]:signal[1]>>hHe3", "he3_cut", "COLZ");
he3->Draw("L SAME");
c1->Draw();
Long64_t n_gate = tree_in->GetEntries("he3_cut");
Long64_t n_true = tree_in->GetEntries("he3_cut && A==3 && Z==2");
std::cout << "Selected entries: " << n_gate << std::endl;
std::cout << "True 3He among selected: " << n_true << std::endl;
std::cout << "Other isotopes among selected: " << n_gate - n_true << std::endl;
Selected entries: 7422 True 3He among selected: 7416 Other isotopes among selected: 6
这次选择只用了测量信号，随后才用模拟的 A、Z 核对成分。区域只覆盖一段条带，边缘也可能漏选或混入其他粒子，因此所选计数不等于全部 ³He 的产额。
4.4 从条带到 PID 投影
第一章用 PID 参数将弯曲条带变成近似水平带，再作一维投影。对本页的示意模型，由 \(\Delta E_{\rm true}=KZ^2A/E_0\) 和 \(E_0=\Delta E_{\rm true}+E_{\rm true}\)，可定义
\[P=\sqrt{\Delta E(\Delta E+E)}\approx\sqrt{KZ^2A}.\]
无测量涨落时，同一核素的 P 在这个模型中是常数；加入测量涨落后形成有宽度的峰。这个表达式与本例模型相配，并不是实际望远镜数据的通用 PID 公式。
tree_in->SetAlias("pid", "sqrt(signal[0]*Etot)");
TH2F *hPidE = new TH2F("hPidE", "PID transformation;E (MeV);P (MeV)",
425, 0, 85, 320, 0, 16);
tree_in->Draw("pid:signal[1]>>hPidE", "", "goff");
TH1D *hAllP = hPidE->ProjectionY("hAllP");
TH1D *hGateP = new TH1D("hGateP", "", 320, 0, 16);
tree_in->Draw("pid>>hGateP", "he3_cut", "goff");
goff 表示先填充、不立即画图。ProjectionY 沿 x 方向累加二维 bin，得到纵轴 P 的分布。cut 仍使用原来的 ΔE、E 坐标，但绘图量已经换成 P。
TCanvas *cPID = new TCanvas("cPID", "PID and projection", 1100, 460);
cPID->Divide(2, 1);
cPID->cd(1);
hPidE->Draw("COLZ");
cPID->cd(2);
hAllP->SetTitle("PID projection;P (MeV);Counts / 0.05 MeV");
hAllP->SetLineColor(kBlack);
hAllP->SetMaximum(1.35 * hAllP->GetMaximum());
hAllP->Draw("HIST");
hGateP->SetLineColor(kRed + 1);
hGateP->Draw("HIST SAME");
TLegend *legendPID = new TLegend(0.55, 0.76, 0.88, 0.89);
legendPID->AddEntry(hAllP, "All events", "l");
legendPID->AddEntry(hGateP, "3He graphical cut", "l");
legendPID->Draw();
cPID->Draw();
左图显示变换后的条带，右图比较全部事例与 ³He cut 的投影。p、d、t、³He、⁴He、⁶He 的中心依次约为 2.83、4.00、4.90、9.80、11.31、13.86 MeV，来自 \(\sqrt{KZ^2A}\)。红色投影落在 ³He 峰附近；氢同位素投影有明显重叠，不能把每个局部起伏都解释成独立核素。
4.5 数组选择与事例选择
signal[0] 每个 entry 只取一个值；直接使用 signal 则展开数组元素，一个 entry 可能贡献多个值。若要判断整个事例，可用 Sum$、Max$ 等在同一 entry 内汇总数组。
std::cout << "Both positive (explicit): "
<< tree_in->GetEntries("signal[0]>0 && signal[1]>0") << std::endl;
std::cout << "Both positive (array): "
<< tree_in->GetEntries("Sum$(signal>0)==2") << std::endl;
std::cout << "At least one above 40 MeV: "
<< tree_in->GetEntries("Max$(signal)>40") << std::endl;
Both positive (explicit): 59977 Both positive (array): 59977 At least one above 40 MeV: 37990
Sum$(signal>0) 统计本事例中正信号的个数，Max$(signal) 取本事例的最大信号，不是在所有 entry 之间求和或取最大值。三层探测器和波形数组采用同样的规则。
Draw 的 selection 实际是填充权重：逻辑条件得到 0 或 1。若直接写 "signal[0]"，含义是用能量作为权重，而不是“选择这个 branch”。本例均用逻辑条件，直方图内容就是所选事例的计数。
5. 在逐事例循环中复用 cut
复杂的逐事例修正、波形积分等计算更适合显式循环。这里保留同一个 ³He 区域，用 IsInside(x, y) 判断当前事例是否在多边形内，并计算所选粒子的总能量谱。它的参数顺序是 x、y，即 E、ΔE。
Int_t eventID_r;
Int_t A_r;
Int_t Z_r;
Float_t E0_r;
Float_t signal_r[2];
tree_in->SetBranchAddress("eventID", &eventID_r);
tree_in->SetBranchAddress("A", &A_r);
tree_in->SetBranchAddress("Z", &Z_r);
tree_in->SetBranchAddress("E0", &E0_r);
tree_in->SetBranchAddress("signal", signal_r);
SetBranchAddress(branch_name, memory_address) 指定读入 branch 数据时写到哪块内存。名称应与文件完全一致，接收变量的类型、数组容量应与 leaf 描述匹配。与写入时相同，标量需用 &，数组名则已提供首元素地址。
后缀 _r 只是用来区分读取 buffer 和前面的写入 buffer，不是 ROOT 的要求。在交互式 C++ notebook 中，随意重复声明同名变量可能报错。
for (Long64_t i = 0; i < 5; ++i) {
Long64_t bytes_read = tree_in->GetEntry(i);
std::cout << eventID_r << " "
<< A_r << " "
<< Z_r << " "
<< signal_r[0] << " "
<< signal_r[1] << " bytes="
<< bytes_read << std::endl;
}
0 2 1 0.29599 65.9363 bytes=24 1 6 2 7.54924 17.6063 bytes=24 2 6 2 3.4103 54.0821 bytes=24 3 2 1 0.369267 53.2156 bytes=24 4 4 2 3.30148 35.0703 bytes=24
GetEntry(i) 将编号为 i 的记录加载到已连接的 buffer，并返回读取字节数，而不是事例对象。返回值非正表示没有读到 entry 数据。编号使用 Long64_t，因为大型 TTree 的记录数可能超过 32 位整数的表示范围。
先用 Draw 得到对照能谱，再建立相同 bin 的空直方图用于循环。这里两幅能谱都只填入 cut 内的事例，不混用不同粒子的信号。
TH1F *hEnergyDraw = new TH1F("hEnergyDraw", "3He gate;#DeltaE + E (MeV);Counts / MeV", 85, 0, 85);
TH1F *hEnergyLoop = new TH1F("hEnergyLoop", "", 85, 0, 85);
tree_in->Draw("Etot>>hEnergyDraw", "he3_cut", "goff");
Long64_t selected_entries = 0;
for (Long64_t i = 0; i < tree_in->GetEntries(); ++i) {
tree_in->GetEntry(i);
Double_t dE = signal_r[0];
Double_t E = signal_r[1];
if (he3->IsInside(E, dE)) {
hEnergyLoop->Fill(dE + E);
++selected_entries;
}
}
std::cout << "Selected by IsInside: " << selected_entries << std::endl;
Selected by IsInside: 7422
这个循环与 Draw(..., "he3_cut") 应选择完全相同的事例。循环中的条件是普通程序语句：C++ 使用 &&、||、!，Python 使用 and、or、not，区别于前面的 ROOT 字符串表达式。
c1->cd();
hEnergyDraw->SetLineColor(kBlack);
hEnergyDraw->SetMaximum(1.35 * hEnergyDraw->GetMaximum());
hEnergyDraw->Draw("HIST");
hEnergyLoop->SetMarkerStyle(20);
hEnergyLoop->SetMarkerSize(0.6);
hEnergyLoop->SetMarkerColor(kRed + 1);
hEnergyLoop->Draw("P SAME");
TLegend *legendE = new TLegend(0.60, 0.73, 0.88, 0.88);
legendE->AddEntry(hEnergyDraw, "TTree::Draw", "l");
legendE->AddEntry(hEnergyLoop, "Loop + IsInside", "p");
legendE->Draw();
c1->Draw();
Double_t max_difference = 0;
for (Int_t b = 0; b <= hEnergyDraw->GetNbinsX() + 1; ++b) {
Double_t difference = std::abs(hEnergyDraw->GetBinContent(b) - hEnergyLoop->GetBinContent(b));
max_difference = std::max(max_difference, difference);
}
std::cout << "Largest bin-content difference: " << max_difference << std::endl;
Largest bin-content difference: 0
两种实现应逐 bin 一致，比较也包含 underflow 和 overflow。能谱的两端受 graphical cut 覆盖范围限制；这里展示的是所选区域内的总能量分布，不是完整的 ³He 入射能谱。
6. 保存 cut、直方图与所选事例
将所选直方图和 cut 一起保存，下次可以复查选择区域。若后续只分析这些事例，CopyTree(cut) 可复制满足条件的完整 entry，包括所有原有 branch。
TFile *fanalysis = new TFile("telescope_analysis_cpp.root", "RECREATE");
fanalysis->cd();
he3->Write();
hPID->Write();
hAllP->Write();
hGateP->Write();
hEnergyDraw->Write();
TTree *selected_tree = tree_in->CopyTree("he3_cut");
selected_tree->Write("he3");
std::cout << "Saved entries: " << selected_tree->GetEntries() << std::endl;
fanalysis->Close();
Saved entries: 7422
fanalysis.cd() 将新 tree 建立在输出文件中；写出的内部名称为 he3。原始输入 TTree 没有改变。he3_cut 则保存了多边形顶点和坐标表达式，重新读取时仍应核对它与待分析数据的 branch、单位和刻度是否一致。
if (fin) fin->Close();
不再使用 TTree 后再关闭输入文件。从 ROOT 文件中读出的对象可能由该文件管理，关闭后继续使用可能访问无效指针；已经写入独立分析文件的直方图则保存在磁盘上。
7. 作业中的定长数组
同样的 branch 写法也可保存探测器能量或波形采样点：
| 数据来源 | leaf 描述 | 一条 entry 的内容 |
|---|---|---|
| 本例 | signal[2]/F | 两层探测器的信号 |
| 作业 4.1 | e[3]/F, A/I, Z/I | D1、D2、D3 的能量沉积及粒子标签 |
| 作业 5.1 | tree wave 中的 adc[250]/D | 一个脉冲的 250 个采样点 |
/F 对应 Float_t，/D 对应 Double_t，/I 对应 Int_t。内存中数组的类型和长度应与文件中的定义一致。
对于波形，GetEntry(i) 一次加载整个脉冲，再通过内层循环处理采样数组，计算基线或电荷积分。下一次 GetEntry 会用下一个脉冲替换 buffer 中的内容。