Armory/aarootscript.C

@@ -46,10 +46,10 @@ void aarootscript(int argument = 0) {
 
     std::cout << "Zero Count: " << zeroCount << std::endl;
 
-    //std::cout << "Making histograms..." << std::endl;
+    std::cout << "Making histograms..." << std::endl;
 
     gErrorIgnoreLevel = 2001;
-    //gROOT->ProcessLine(".x histcomp.C");
+    gROOT->ProcessLine(".x histcomp.C");
 
     std::cout << "=========================================\n";
 
@@ -68,59 +68,7 @@ void aarootscript(int argument = 0) {
             std::cout << "Tree2 not found!" << std::endl;
         }
     }
-
-    std::cout << "Creating Tb vs dEb plot..." << std::endl;
-
-    // Readers for both trees
-    TTreeReader r1(tree1);
-    TTreeReader r2(tree2);
-
-    TTreeReaderValue<double> Tb_val(r1, "Tb");
-    TTreeReaderValue<double> dEb_val(r2, "dEb");
-
-    std::vector<double> x; // Tb (tree1)
-    std::vector<double> y; // dEb (tree2)
-
-
-    // Loop over both trees simultaneously
-    while (r1.Next() && r2.Next()) {
-        x.push_back(*Tb_val);
-        y.push_back(*dEb_val);
-    }
-    std::cout << "x length: " << x.size() << ", y length: " << y.size() << std::endl;
-
-    #include <fstream>
-
-    std::ofstream outfile("Tb_dEb_data.txt");
-
-    if (!outfile.is_open()) {
-        std::cerr << "Error: Could not open output file!" << std::endl;
-        return;
-    }
-
-    for (size_t i = 0; i < x.size(); i++) {
-        outfile << x[i] << " " << y[i] << "\n";
-    }
-
-    outfile.close();
-
-    std::cout << "Data written to Tb_dEb_data.txt" << std::endl;
-    /*
-    // Create graph
-    TGraph *gr = new TGraph(x.size(), &x[0], &y[0]);
-    gr->SetTitle("Tb (tree1) vs dEb (tree2);Tb;dEb");
-    gr->SetMarkerStyle(20);
-
-    // Draw
-    TCanvas *c1 = new TCanvas("c1", "Tb vs dEb", 800, 600);
-    gr->Draw("AP");
-    c1->Update();
-    c1->SaveAs("Tb_vs_dEb.png");
-    std::cout << "Plot saved as Tb_vs_dEb.pdf" << std::endl;
     std::cout << "\n\n\n";
 
-    delete c1;
-    delete gr;*/
-
 }
 

Armory/anasenMS.cpp

@@ -35,7 +35,7 @@ int main(int argc, char **argv){
   // load energy loss tables (assume units: E in MeV, dE/dx in MeV/(mg/cm²), density in mg/cm³)
   TGraph* elossLight = LoadELoss("../ELoss/Eloss_HeAlpha"); // for light particle (alpha)
   TGraph* elossHeavy = LoadELoss("../ELoss/Eloss_p");     // for heavy particle (proton)
-  double density = (2.1525e-7) * 400; // example for aluminum target, adjust as needed (400 torr is 0.0000861)
+  double density = 0.0000861; // example for aluminum target, adjust as needed
   auto c1 = new TCanvas("c1", "Graph Example", 800, 600);
   auto g = elossLight;
   g->SetTitle("Energy Loss Table;Kinetic Energy (MeV);dE/dx (MeV/(mg/cm^{2}))");
@@ -304,10 +304,11 @@ int main(int argc, char **argv){
         printf("Event %d: Ekin before loss = %f MeV, distance = %f cm\n", i, Tb, dl);
       }
       double tb_temp = Tb;
+      //double dE_light = dedxLight * dl * density / 1000.0; // adjust for units (example scaling)
       double dx = 0;
       double counter = 0;
       while(dx < dl){
-        double step = 0.1; // cm, step size for energy loss calculation
+        double step = 0.01; // cm, step size for energy loss calculation
         if(dx + step > dl) step = dl - dx; // adjust last step to end at hit position
         double EkinStep = Tb; // kinetic energy at current step
         double dedxStep = elossLight->Eval(EkinStep); // dE/dx at current energy