Hi all,
My code is following
#include <string>
#include <sstream>
#include <fstream>
#include <unistd.h>
#include <stdlib.h>
#include <algorithm>
#include "struct.h"
#include "ReadGDML.h"
#include "ReadEDepSim.h"
#include "ReadOptical.h"
#include "ReadDigit.h"
#include "CameraImage.h"
#include "FastFitter.h"
#include "LensDisplayer.h"
#include "LineIntersector.h"
#include "StructDefinition.h"
#include "TrackMatcher.h"
#include "Volume3DReconstructer.h"
#include "EnergyEstimator.h"
#include "TFile.h"
#include "TTree.h"
#include "TMath.h"
#include "TApplication.h"
#include "TF1.h"
#include "TH2I.h"
#include "TCanvas.h"
#include "TVector3.h"
#include "TH1.h"
#include "TRandom.h"
#include "TSpectrum.h"
#include "TChain.h"
#include "TVirtualFitter.h"
#include "Riostream.h"
int npeaks=10 ;
int final_nfound;
int final_nfound1;
int A;
int B;
int t;
float taufast=6e-9;
float tauslow=1.6e-6;
void usage() {
std::cout << "Input parameters needed: " << std::endl;
std::cout << "-i 'root file with spill data '" << std::endl;
std::cout << std::endl;
}
int main(int argc, char **argv) {
std::string inFilename="";
std::string outFilename="";
int iflag=0, oflag=0; //mandatory
int opt;
outFilename = "interaction.root";
while( (opt = getopt(argc, argv, "i:o")) != -1) {
switch (opt) {
case 'o':
if(optarg) outFilename = optarg;
oflag = 1;
break;
case 'i':
if(optarg) inFilename = optarg;
iflag = 1;
break;
default:
usage();
return 1;
}
}
if (!iflag) {
std::cout << "ERROR: missing REQUIRED parameters!" << std::endl;
usage();
return 1;
}
std::cout << "Sensor file: " << inFilename << std::endl;
std::cout << "Output file: " << outFilename << std::endl;
std::cout << "*****************************************************" << std::endl;
std::cout << "Reading data" << std::endl;
TFile *m_f = new TFile(inFilename.c_str(), "READ");
int k;
TChain *tree = new TChain();
for (int k=1;k<39;k++){
tree->Add(Form("./%s/CAM_%d",inFilename.c_str(),k));
}
if (!tree) {
std::cout << "ERROR: missing the event trees!" << std::endl;
}
int idEvent;
std::vector<double> *energy=0 ;
std::vector<double> *time=0 ;
tree->SetBranchAddress("time", &time);
tree->SetBranchAddress("idEvent",&idEvent);
tree->SetBranchAddress("energy",&energy);
const int n=9;
TH1F *peakno[n];
TH1F *intno[n];
TH1F *inttotal= new TH1F ("inttotal","inttotal",100,0,100);
for (int i=0;i<n;i++) {
peakno[i] = new TH1F(Form("peakno%d",i),Form("peakno%d",i),20000,0,20000);
intno[i] = new TH1F(Form("intno%d",i),Form("intno%d",i),100,0,100);
}
int N = tree->GetEntries();
for (int i=0; i<N; i++) {
tree->GetEntry(i);
ULong_t nsig = time->size();
double *sig = time->data();
double *esig = energy->data();
if (idEvent>=0 && idEvent<=8) {
for (ULong_t j=0; j<nsig; j++) {
if (esig[j] < 6.0e-6 ) continue;
peakno[idEvent]->Fill(sig[j]);
}
}
}
for (int i=0;i<n;i++) {
TSpectrum *s = new TSpectrum(npeaks);
int nfound = s->Search(peakno[i],0.0,"",0.1);
intno[i]->Fill(nfound);
std::cout << "idEvent = " << i << "no of peaks= " << nfound << std::endl;
inttotal->Add(intno[i]);
}
TFile *outputFile = new TFile(outFilename.c_str(),"RECREATE");
std::cout << "Output file written to " << outFilename << std::endl;
for (int i=0;i<n;i++) {
peakno[i]->Write();
intno[i]->Write();
}
inttotal->Write();
outputFile->Close();
}
After running the code i generated few histogram with no of peaks. The output is like,
Reading data
Warning in <TSpectrum::SearchHighRes>: Peak buffer full
Info in <TCanvas::MakeDefCanvas>: created default TCanvas with name c1
idEvent = 0no of peaks= 10
idEvent = 1no of peaks= 2
idEvent = 2no of peaks= 3
idEvent = 3no of peaks= 1
idEvent = 4no of peaks= 1
idEvent = 5no of peaks= 5
idEvent = 6no of peaks= 2
idEvent = 7no of peaks= 3
idEvent = 8no of peaks= 1
Output file written to interaction.root
Now in the histogram,
peakno0.pdf (80.3 KB)
there is no actual peak…There is no counts… I want to put a thresold into the yaxis count in the code so that this histogram will not be generated. Is it possible to give yaxis count before generating the histogram into this code ?
If we see another example of histogram,
peakno1.pdf (111.7 KB)
this shows two nice peaks.
Thank you.