Sorry for the late reply and thanks for your patience. I can not even compile my macro, because when I try to do
.x name_of_macro.C+ I get
root [0] .x trackingchi2.C+
Info in <TUnixSystem::ACLiC>: creating shared library /home/matteo/mu_beam bumps/./trackingchi2_C.so
start processing
*** Break *** segmentation violation
Similarly to when I try to run the macro from root. If I write a c++ program and I compile it with c++ compiler, it compiles but I get the same “segmentation violation” error when I run it. While if I compile it with gcc I get:
gcc -o macro macro.cpp `root-config --cflags --glibs`
/usr/bin/ld: /tmp/ccQauRBx.o: undefined reference to symbol '_ZSt20__throw_length_errorPKc@@GLIBCXX_3.4'
/usr/lib/libstdc++.so.6: error adding symbols: DSO missing from command line
collect2: error: ld returned 1 exit status
The content of the macro is the following:
#include "TFile.h"
#include "TTree.h"
#include "TGraph.h"
#include "TGraphErrors.h"
#include <vector>
#include <algorithm>
#include "TH1F.h"
#include "TH2F.h"
#include "TCanvas.h"
#include "TF1.h"
#include <iostream>
#include <fstream>
#include <string>
#include <cmath>
#include "TMath.h"
using namespace std;
int main(int argc, char* argv[]) {
ofstream myfile;
const char* label = "trackingchi2";
TString fileName = "trackingTB_run5367_68_70_71_72_73_74new.root";
myfile.open("trackingchi2.txt");
TFile fhisto( TString(label) + ".root","RECREATE");
TGraph* g_theta = new TGraph();
g_theta->SetTitle("#theta_{e} vs #theta_{mu}");
g_theta->SetMarkerStyle(20);
int run;
int eventID;
int SingleTrack;
int signal;
double chi2_x, chi2_y;
std::vector<double> *hits0x, *hits0y, *hits1x, *hits1y, *hits2x, *hits2y, *hits3x, *hits3y, *hits4x, *hits4y;
std::vector<double> *thetaIn_x, *thetaIn_y, *thetaOut_x, *thetaOut_y;
TFile * file = new TFile(fileName);
TTree * tree = (TTree*)file->Get("simpleEvent");
Int_t nentries = (Int_t)tree->GetEntries();
TBranch * brun = tree->GetBranch("run");
TBranch * beventID = tree->GetBranch("eventID");
TBranch * bSingleTrack = tree->GetBranch("SingleTrack");
TBranch * bsignal = tree->GetBranch("signal");
TBranch * bHits0x = tree->GetBranch("Hits0x");
TBranch * bHits0y = tree->GetBranch("Hits0y");
TBranch * bHits1x = tree->GetBranch("Hits1x");
TBranch * bHits1y = tree->GetBranch("Hits1y");
TBranch * bHits2x = tree->GetBranch("Hits2x");
TBranch * bHits2y = tree->GetBranch("Hits2y");
TBranch * bHits3x = tree->GetBranch("Hits3x");
TBranch * bHits3y = tree->GetBranch("Hits3y");
TBranch * bHits4x = tree->GetBranch("Hits4x");
TBranch * bHits4y = tree->GetBranch("Hits4y");
TBranch * bthetaIn_x = tree->GetBranch("thetaIn_x");
TBranch * bthetaIn_y = tree->GetBranch("thetaIn_y");
TBranch * bthetaOut_x = tree->GetBranch("thetaOut_x");
TBranch * bthetaOut_y = tree->GetBranch("thetaOut_y");
TBranch * bchi2_x = tree->GetBranch("chi2_x");
TBranch * bchi2_y = tree->GetBranch("chi2_y");
brun->SetAddress(&run);
beventID->SetAddress(&eventID);
bSingleTrack->SetAddress(&SingleTrack);
bsignal->SetAddress(&signal);
bHits0x->SetAddress(&hits0x);
bHits0y->SetAddress(&hits0y);
bHits1x->SetAddress(&hits1x);
bHits1y->SetAddress(&hits1y);
bHits2x->SetAddress(&hits2x);
bHits2y->SetAddress(&hits2y);
bHits3x->SetAddress(&hits3x);
bHits3y->SetAddress(&hits3y);
bHits4x->SetAddress(&hits4x);
bHits4y->SetAddress(&hits4y);
bthetaIn_x->SetAddress(&thetaIn_x);
bthetaIn_y->SetAddress(&thetaIn_y);
bthetaOut_x->SetAddress(&thetaOut_x);
bthetaOut_y->SetAddress(&thetaOut_y);
bchi2_x->SetAddress(&chi2_x);
bchi2_y->SetAddress(&chi2_y);
// lines for geom configuration description
double z0= 0.0 ;
double z1= 9390.0 ;
double Dz1= z1 - z0 ;
// measured z of entering face of the target
double zT= 9709.0 ;
double zT_mid= 9709.0 + 4 ;
// z position from logbook if the 3 downstream planes
double z2= 9788.0 ;
double z3= 10290.0 ;
double z4= 10698.0 ;
const double E_mu_160 = 160;
const double E_mu_176 = 176;
const double E_mu_144 = 144;
const double m_mu = 0.1056583745;
const double m_e = 0.000511;
const double rsq_160 = (pow(E_mu_160,2)-pow(m_mu,2))/pow(E_mu_160+m_e,2);
TF1* elasticity_160 = new TF1("elasticity_160", "asin(sin(x)*([0]+[0]*[1]*pow(cos(x),2))/([2]-2*[0]*[1]*pow(cos(x),2)-[2]*[1]*pow(cos(x),2)))", 0., 0.05);
elasticity_160->SetParameter(0, m_e);
elasticity_160->SetParameter(1, rsq_160);
elasticity_160->SetParameter(2, E_mu_160);
elasticity_160->SetLineColor(kRed);
double DzT4= z4 - zT_mid ;
double Dz34 = z4- z3 ;
double Dz24 = z4- z2 ;
double DzT3 = DzT4-Dz34;
double DzT2 = DzT4-Dz24;
double zT_mid_0 = zT_mid-z1;
double z2_0 = z2-z1;
double z3_0 = z3-z1;
double z4_0 = z4-z1;
std::vector<double> chi2_2x, chi2_2y, chi2_3x, chi2_3y, best_2x, best_2y, best_3x, best_3y, best_4x, best_4y;
std::vector<double> chi2_4x, chi2_4y, chi2_best, real_2x, real_2y, real_3x, real_3y, x_4, y_4, theta_x, theta_y, angle_x, angle_y;
double best_theta_e_x, best_theta_mu_x, best_theta_e_y, best_theta_mu_y, min_chi2;
int graph_index = 1;
std::cout << "start processing" << std::endl;
double mythin_x = 999.;
double mythin_y = 999.;
double xthout_4;// = 999.;
double ythout_4;// = 999.;
double xthout_3 = 999.;
double ythout_3 = 999.;
//double best_2x, best_2y, best_3x, best_3y, best_chi2x, best_chi2y;
double scarto_2x;// = 90;
double scarto_3x;// = 90.;
double scarto_2y;// = 90.;
double scarto_3y;// = 90.;
int minPos = 0;
int minPos2 = 0;
int maxPos = 0;
for (Int_t i = 0; i < nentries; i++) {
Long64_t tentry = tree->LoadTree(i);
tree->GetEntry(i);
// selects events with only 1 hit in pl0 and pl1 (incoming track)
if ((hits0x->size() * hits0y->size() *
hits1x->size() * hits1y->size()) == 1) {
if(hits4x->size()==4 && hits4y->size()==4 && hits3x->size()==4 && hits3y->size()==4 && hits2x->size()==4 && hits2y->size()==4){ // && hits4y->size()>=2 && hits4y->size()<=3) {
double dx = fabs(hits4x->at(1)-hits4x->at(0));
double dy = fabs(hits4y->at(1)-hits4y->at(0));
double x0pl0 = hits0x->at(0) ;
double x1pl1 = hits1x->at(0) ;
double y0pl0 = hits0y->at(0) ;
double y1pl1 = hits1y->at(0) ;
mythin_x = ((x1pl1 - x0pl0 ) / Dz1 ) ;
mythin_y = ((y1pl1 - y0pl0 ) / Dz1 ) ;
// expected x,y coord in the target
Double_t xin_target = x1pl1 + mythin_x * (zT_mid-z1);
Double_t yin_target = y1pl1 + mythin_y * (zT_mid-z1);
double thetaIn = (sqrt(pow(x1pl1-x0pl0,2)+pow(y1pl1-y0pl0,2))/Dz1);
for (unsigned j = 0; j < hits2x->size(); j++) {
for (unsigned l = 0; l < hits3x->size(); l++) {
for (unsigned m = 0; m < hits4x->size(); m++) {
TGraphErrors* g_x = new TGraphErrors();
TF1* fitf_x = new TF1("fitf_x", "pol1", z2, z4);
g_x->SetPoint(0, z2, hits2x->at(j));
g_x->SetPoint(1, z3, hits3x->at(l));
g_x->SetPoint(2, z4, hits4x->at(m));
g_x->SetPointError(0, 0, 0.007);
g_x->SetPointError(1, 0, 0.007);
g_x->SetPointError(2, 0, 0.035);
g_x->Fit("fitf_x", "Q");
double dx2 = pow((fitf_x->Eval(z2)-hits2x->at(j))/(0.007),2);
double dx3 = pow((fitf_x->Eval(z3)-hits3x->at(l))/(0.007), 2);
double dx4 = pow((fitf_x->Eval(z4)-hits4x->at(m))/(0.035), 2);
chi2_3x.push_back(dx2+dx3+dx4);
theta_x.push_back((fitf_x->GetParameter(1)));
best_2x.push_back(hits2x->at(j));
best_3x.push_back(hits3x->at(l));
best_4x.push_back(hits4x->at(m));
g_x->Delete();
fitf_x->Delete();
}
}
}
for (unsigned j = 0; j < hits2y->size(); j++) {
for (unsigned l = 0; l < hits3y->size(); l++) {
for (unsigned m = 0; m < hits4y->size(); m++) {
TGraphErrors* g_y = new TGraphErrors();
TF1* fitf_y = new TF1("fitf_y", "pol1", z2, z4);
g_y->SetPoint(0, z2, hits2y->at(j));
g_y->SetPoint(1, z3, hits3y->at(l));
g_y->SetPoint(2, z4, hits4y->at(m));
g_y->SetPointError(0, 0, 0.007);
g_y->SetPointError(1, 0, 0.007);
g_y->SetPointError(2, 0, 0.035);
g_y->Fit("fitf_y", "Q");
double dy2 = pow((fitf_y->Eval(z2)-hits2y->at(j))/(0.007),2);
double dy3 = pow((fitf_y->Eval(z3)-hits3y->at(l))/(0.007), 2);
double dy4 = pow((fitf_y->Eval(z4)-hits4y->at(m))/(0.035), 2);
chi2_3y.push_back(dy2+dy3+dy4);
theta_y.push_back((fitf_y->GetParameter(1)));
best_2y.push_back(hits2y->at(j));
best_3y.push_back(hits3y->at(l));
best_4y.push_back(hits4y->at(m));
g_y->Delete();
fitf_y->Delete();
}
}
}
for (unsigned j = 0; j < chi2_3x.size(); j++) {
chi2_best.push_back((chi2_3x[j]+chi2_3y[j]));
}
min_chi2 = 1e5;
for (unsigned a = 0; a < chi2_best.size(); a++) {
if (chi2_best[a] < min_chi2) {
min_chi2 = chi2_best[a];
minPos = a;
}
}
min_chi2 = 1e5;
for (unsigned a = 0; a < chi2_best.size(); a++) {
if (chi2_best[a] == chi2_best[minPos]) continue;
if (chi2_best[a] < min_chi2) {
min_chi2 = chi2_best[a];
minPos2 = a;
}
}
best_theta_mu_x = theta_x[minPos];
best_theta_e_x = theta_x[minPos2];
best_theta_mu_y = theta_y[minPos];
best_theta_e_y = theta_y[minPos2];
double thetaE = acos((1+mythin_x*best_theta_e_x+mythin_y*best_theta_e_y)/sqrt((1+mythin_x*mythin_x+mythin_y*mythin_y)*(1+best_theta_e_x*best_theta_e_x+best_theta_e_y*best_theta_e_y)));
double thetaMU = acos((1+mythin_x*best_theta_mu_x+mythin_y*best_theta_mu_y)/sqrt((1+mythin_x*mythin_x+mythin_y*mythin_y)*(1+best_theta_mu_x*best_theta_mu_x+best_theta_mu_y*best_theta_mu_y)));
g_theta->SetPoint(graph_index, thetaE, thetaMU);
++graph_index;
chi2_3x.clear();
chi2_3y.clear();
theta_x.clear();
theta_y.clear();
chi2_best.clear();
best_2x.clear();
best_3x.clear();
best_4x.clear();
best_2y.clear();
best_3y.clear();
best_4y.clear();
}// end over selection > 1 hit in pl4
}// end loop over selection = 1 hit in pl0 pl1 only
if(i%100000 == 0) std::cout << "...... processing" << std::endl;
} // end loop over events
std::cout << "processing finished" << std::endl;
TCanvas* c2 = new TCanvas();
g_theta->Draw("A*");
g_theta->GetYaxis()->SetRangeUser(0., 0.005);
g_theta->GetXaxis()->SetRangeUser(0., 0.05);
g_theta->SetMarkerSize(0.7);
elasticity_160->SetNpx(600);
elasticity_160->Draw("SAME");
c2->Modified();
c2->Update();
c2->Print("Tchi2thetaemu.root", "root");
fhisto.Write();
myfile.close();
return 0;
}
Unfortunately the upload of the TTree I am trying to read systematically fails…
Regards,
Matteo
).