前述分析程序的做法有两个主要的缺点：
避免上述问题的更好的方法是采用分析代码类(class ana)继承MakeClass生成的类(class tracking)的方式。通过这样的方式可不改动原来的类，将用户分析代码与原来的代码分开，互不干扰。
主目录: ./tracking
./tracking
- main.cpp - 主程序
- Makefile -编译方法
./tracking/include/
- tracking.h - MakeClass 生成的头文件
- ana.h - 用户分析头文件
./tracking/src/
- tracking.C -MakeClass生成的源文件
- ana.cpp 用户分析源文件#include <iostream>
#include <TFile.h>
#include <TTree.h>
#include <TString.h>
#include "ana.h" //用户代码头文件
using namespace std;
int main(int argc, char** argv)
{
if(argc !=2) {
cout<<"Usage: ./"<<argv[0]<<" run_number "<<endl;
return -1;
}
int run_number = atoi(argv[1]);
TString InputPath, OutputPath, infile, outfile;
InputPath = "./";
OutputPath = "./";
infile.Form("%sf8ppac%03d.root", InputPath.Data(), run_number);
outfile.Form("%sout%03d.root", OutputPath.Data(), run_number);
//input
TFile *ipf = new TFile(infile);
if(!ipf->IsOpen()) {
cout<<"Cannot open input file: "<<infile<<endl;
return -1;
}
TTree *ipt = (TTree*)ipf->Get("tree");
//output
TFile *opf = new TFile(outfile,"RECREATE");
TTree *opt = new TTree("tree","ppac tracking");
//
ana *tk = new ana(run_number,ipt);
tk->Analysis(opt); //等价于原Loop函数
//
ipf->Close();
opf->Close();
return 1;
}
#ifndef ana_h
#define ana_h
#include <TH2.h>
#include <iostream>
#include "tracking.h" //包含基类头文件
using namespace std;
class ana : public tracking //从tracking类中继承其成员变量和成员函数
{
public:
int run; //run_number
Double_t xx[3],xz[3],yy[3],yz[3];//1A,2A,3
Double_t p2bx,p2by;//2B x-y
Double_t dx[3],dy[3];
Double_t tx,ty;//target pos
Double_t c2nx,c2ny;//chi2/ndf for xfit,yfit
TH2D *htf8x,*htf8y,*htar;
//构造函数：将传入的tree传递给tracking基类，完成初始化
ana(int run_number, TTree* tree): tracking(tree), run(run_number) {} //run=run_number
virtual ~ana() {};
virtual void Analysis(TTree *tree);//分析函数，作用等价于原Loop函数
virtual void SetOutBranch(TTree *tree);
virtual void TrackInit();
virtual void SetTrace(TH2D *h,Double_t k,Double_t b,Int_t min,Int_t max);
};
#endif
#include <TH2.h>
#include <TStyle.h>
#include <TCanvas.h>
#include <TF1.h>
#include <TFitResult.h>
#include <TGraph.h>
#include "ana.h"
using namespace std;
void ana::SetOutBranch(TTree *tree)
{
//measured pos
tree->Branch("xx",&xx,"xx[3]/D");//1A,2A,3
tree->Branch("xz",&xz,"xz[3]/D");
tree->Branch("yy",&yy,"yy[3]/D");
tree->Branch("yz",&yz,"yz[3]/D");
//difference between measured and calculated -for pos resolution.
tree->Branch("dx",&dx,"dx[3]/D");
tree->Branch("dy",&dy,"dy[3]/D");
//for efficiency calculation -2B
tree->Branch("p2bx",&p2bx,"p2bx/D");
tree->Branch("p2by",&p2by,"p2by/D");
//target x-y
tree->Branch("tx",&tx,"tx/D");
tree->Branch("ty",&ty,"ty/D");
//ch2/ndf for linear fitting.
tree->Branch("c2nx",&c2nx,"c2nx/D");
tree->Branch("c2ny",&c2ny,"c2ny/D");
}
void ana::TrackInit()
{
tx=-999;
ty=-999;
//1A
xx[0]=PPACF8[0][0];
yy[0]=PPACF8[0][1];
xz[0]=PPACF8[0][2];
yz[0]=PPACF8[0][3];
//2A
xx[1]=PPACF8[2][0];
yy[1]=PPACF8[2][1];
xz[1]=PPACF8[2][2];
yz[1]=PPACF8[2][3];
//3
xx[2]=PPACF8[4][0];
yy[2]=PPACF8[4][1];
xz[2]=PPACF8[4][2];
yz[2]=PPACF8[4][3];
//2B
p2bx=PPACF8[3][0];
p2by=PPACF8[3][1];
}
void ana::SetTrace(TH2D *h,Double_t k,Double_t b,Int_t min,Int_t max){
if(h==0) return;
if(min>=max) return;
for(int i=min;i<max;i++){
h->Fill(i,(Int_t)(i*k+b));
}
}
void ana::Analysis(TTree *tree)
{
TH2D *htf8x=new TH2D("htf8x","x trace by ppac",2200,-2000,200,300,-150,150);
TH2D *htf8y=new TH2D("htf8y","y trace by ppac",2200,-2000,200,300,-150,150);
TH2D *htar=new TH2D("htar","distribution on target",100,-50,50,100,-50,50);
SetOutBranch(tree);
if (fChain == 0) return;
Long64_t nentries = fChain->GetEntriesFast();
Long64_t nbytes = 0, nb = 0;
for (Long64_t jentry=0; jentry<nentries;jentry++) {
Long64_t ientry = LoadTree(jentry);
if (ientry < 0) break;
nb = fChain->GetEntry(jentry); nbytes += nb;
TrackInit();
bool b1a=xx[0]>-999 && yy[0]>-999;
bool b2a=xx[1]>-999 && yy[1]>-999;
bool b3=xx[2]>-999 && yy[2]>-999;
if(!b1a || !b2a || !b3) continue;
//fit x-z trajectory
TFitResultPtr r;
TGraph *grx=new TGraph(3,xz,xx);
TF1 *fx=new TF1("fx","pol1",-2000,0);
r=grx->Fit(fx,"SQ");
tx=fx->Eval(0);
SetTrace(htf8x,fx->GetParameter(1),fx->GetParameter(0),-1800,0);
for(int i=0;i<3;i++) dx[i]=xx[i]-fx->Eval(xz[i]);
c2nx=r->Chi2()/r->Ndf();
delete grx;
delete fx;
//fit y-z trajectory
TGraph *gry=new TGraph(3,yz,yy);
TF1 *fy=new TF1("fy","pol1",-2000,0);
r=gry->Fit(fy,"SQ");
ty=fy->Eval(0);
SetTrace(htf8y,fy->GetParameter(1),fy->GetParameter(0),-1800,0);
for(int i=0;i<3;i++) dy[i]=yy[i]-fy->Eval(yz[i]);
c2ny=r->Chi2()/r->Ndf();
delete gry;
delete fy;
htar->Fill(tx,ty);
tree->Fill();
if(jentry%10000==0) cout<<"processing "<<jentry<<endl;
}
htf8x->Write();
htf8y->Write();
htar->Write();
tree->Write();
}
!jupyter nbconvert exercise2.4 --to html