部門/ROOT/キャリブレーションしてROOTファイルをまるまる書き換える - 阪大理学部科学計算・開発サークル -Circle for Science Caluculation and Development-

イベントのごとにデータをすべて書き換えるのは、以前かなり時間がかかると言ったがそれなりの現時点で平均的なデスクトップPCくらいのスペックを持ったものなら、数GBのものを処理するのにも数分程度で終わることが分かった。一度キャリブレーションしておけば、のちのち解析でカットを掛けたりするときにエネルギーや時間で直接条件が指定できるので、便利だろうということで、TTreeとTFIleを勉強してつくってみた。K本研のDAQシステムで得られたデータを、K本研に伝わるプログラムを使って変換したROOTファイルに対するマクロなので、あまり参考にならないかもしれない。

#include<fstream>
#include<TH1.h>
#include<TH2.h>
#include<TFile.h>
#include<TTree.h>
#include<TCanvas.h>
#include<parameter20111003.h>
#include<math.h>
//#include<iostream>
 
#include<iostream> using namespace std;
 
 
void calb_rewrite(){
  /*
ADCとTDCのヒストグラムの横軸をエネルギーの時間に変えてデータを全部作りなおすマクロ。
   */
 
 
  //  gStyle->SetOptFit(0001);
  gStyle->SetOptFit(0000);
 
  char inname[100],outname[100];
  sprintf(inname,"22na-10.03.001-002.root");//入力ファイル名
  //  sprintf(inname,"caf2-20111005-21.root");//入力ファイル名
  sprintf(outname,"caf2_test.root");//出力ファイル名
 
  TFile *f1=new TFile(inname);//opening root file
 
 
  //TTree *tree =(TTree*)f1->Get("h40");
  TTree *tree;
  tree=h40;
 
  Int_t nevent;//イベント数
  nevent=tree->GetEntries();
 
  Float_t Naiadc[6];
  tree->SetBranchAddress("Naiadc",&Naiadc);
 
  Float_t Naitdc[6];
  tree->SetBranchAddress("Naitdc",&Naitdc);
 
  Double_t realADC[7],realTDC[7];
  char name[10],title[20];
  TObjArray Hlist(0);
  TH1F* hhh;
 
 
  //  nevent=10000;  //for testing  
 
 
  TFile f(outname,"recreate");
 
//making ADC Histgram
  for (Int_t id =1; id<7; id++){
    sprintf(name,"h10%d",id);
    sprintf(title,"Calibrated ADC ch%d",id);
    hhh = new TH1F(name,title,2000,0,4000);
    Hlist.Add(hhh);
    for(Int_t nev=1;nev<nevent;nev++){
      if (nev%1000000 == 0)
	printf("ADC:id =%d, %dth loop, and nevent =%d\n",id,nev,nevent);
      tree->GetEntry(nev);
 
      realADC[id]=Naiadc[id-1]+para1[id]+para2[id];
      hhh->Fill(realADC[id]);
    }  
  }
 
 
  //making TDC Histgram
  for (Int_t id =1; id<7; id++){
    sprintf(name,"h20%d",id);
    sprintf(title,"Calibrated TDC ch%d",id);
    hhh = new TH1F(name,title,2000,0,4000);
    Hlist.Add(hhh);
    for(Int_t nev=1;nev<nevent;nev++){
      if (nev%1000000 == 0)
	printf("TDC:id =%d, %dth loop, and nevent =%d\n",id,nev,nevent);
      tree->GetEntry(nev);
      realTDC[id]=Naitdc[id-1]+para3[id]+para4[id];
      hhh->Fill(realTDC[id]);
    }  
  }
 
  //making Coincidence Tree
  TTree *tree2= new TTree("h402","Coincidence");
 
  Float_t Naiadc1[6],Naitdc1[6];
   tree2->Branch("Naiadc1",&Naiadc1,"Naiadc1[6]/F");
   tree2->Branch("Naitdc1",&Naitdc1,"Naitdc1[6]/F");
 
   for(Int_t nev=1;nev<nevent;nev++){
     if (nev%1000000 == 0)
       printf(" %dth loop, and nevent =%d\n",nev,nevent);
     tree->GetEntry(nev);
     for (Int_t id=0;id<7;id++){
       Naiadc1[id]=Naiadc[id]*para1[id+1]+para2[id+1];
       Naitdc1[id]=Naitdc[id]*para3[id+1]+para4[id+1];
     }
     tree2->Fill();
   }  
   tree2.Write();
 
   Hlist.Add(tree2);
   Hlist->Write();
   f.Close();
}
 

なお、読み込んでいるparameter.hというファイルは、

Double_t para1[7], para2[7],para3[7],para4[7];
Double_t para12,para13,para14,para15,para16, para23,para24,para25,para26,para34,para35,para36,para45,para46,para56;
 
//Energy = Naiadc[ch]*para1[ch]+para2[ch]
para1[1]=0.807352;
para2[1]=-180.236937;
para1[2]=1.008703;
para2[2]=-335.023608;
para1[3]=0.814718;
para2[3]=-221.704268;
para1[4]=1.148811;
para2[4]=-96.457329;
para1[5]=1.271699;
para2[5]=-432.295675;
para1[6]=0.855447;
para2[6]=-71.29883;
//Time=Naitdc[ch]*para3[ch]+para[4]
para3[1]=0.88;
para4[1]=342.85;
para3[2]=1.76;
para4[2]=-83.28;
para3[3]=1.86;
para4[3]=-124.22;
para3[4]=1.66;
para4[4]=-19.17;
para3[5]=1.77;
para4[5]=-88.49;
para3[6]=1.77;
para4[6]=-80.12;
 
para12=0;
para13=0;
para14=0;
para15=0;
para16=0;
 
para23=0;
para24=0;
para25=0;
para26=0;
 
para34=0;
para35=0;
para36=0;
 
para45=0;
para46=0;
 
para56=0;
 

という感じになっている。para12とかは今は使ってないけどTDCの時間差のスペクトルをつくるのにそのうち使う予定。