Comp Sci Drawing a graph using C++ / ROOT (CERN)

AI Thread Summary
The discussion revolves around difficulties in plotting a graph using C++ and the ROOT framework from CERN. The user is trying to visualize energy deposition values associated with specific x-intervals but is unsure how to implement this using their existing code. Key issues include inconsistent x-interval definitions that could lead to double counting of energy values and confusion regarding the appropriate graphing utility (TGraph vs. TGraph2D). Suggestions emphasize the need to refer to the ROOT documentation for guidance on graph creation and ensuring consistent interval definitions. The user seeks further assistance to refine their approach and successfully generate the desired graph.
lavster
Messages
213
Reaction score
0

Homework Statement



Im having difficulty drawing a specific graph. I am not entirely sure how to explain my question/ the graph i need to draw so i think the best way is to provide a bit of my code and explain using it the graph using it...

the code below gives a read out of the energy deposited (Phantom_Ed) for a particular x interval (Phantom_xpre). These read-outs are given a label - E1 for between -100 and -90mm, E2 for -90 and -80mm etc. What id like to do is produce a graph of the position (Phantom_xpre) with the corresponding value of Energy deposited (E1 for first interval etc) plotted on the y axis. However I am not sure how I am to go about doing this...

Homework Equations



my code at the moment (with the bit before defining the tree etc missing) is:

Long64_t nentries = T->GetEntries();// so can determine how many entries for loop
float E1=0; // defining variables for loop and setting to zero
float E2=0;
float E3=0;
float E4=0;
float E11=0;
float E12=0;
float E13=0;
float E14=0;

for(Long64_t k=0;k<nentries;k++)//original loop - gets the entry
{

T->GetEntry(k);

for( Long64_t l=0; l < Phantom_Nhits; l++ ) //steps loop
{
if((Phantom_xpre[l]>=-100) && (Phantom_xpre[l]<-90)) //if statement defining the region of x-axis of interest
{
E1=Phantom_Ed[l]+E1; //adding up the total energy depsoited in this region
}
if((Phantom_xpre[l]>=-90) && (Phantom_xpre[l]<=-80))
{
E11=Phantom_Ed[l]+E11;
}
if((Phantom_xpre[l]>=-80) && (Phantom_xpre[l]<=-70))
{
E12=Phantom_Ed[l]+E12;
}
if((Phantom_xpre[l]>=-70) && (Phantom_xpre[l]<=-60))
{
E13=Phantom_Ed[l]+E13;
}
if((Phantom_xpre[l]>=-60) && (Phantom_xpre[l]<=-50))
{
E14=Phantom_Ed[l]+E14;
}
if((Phantom_xpre[l]>=-50) && (Phantom_xpre[l]<0))
{
E2=Phantom_Ed[l]+E2;
}
if((Phantom_xpre[l]>=0) && (Phantom_xpre[l]<50))
{
E3=Phantom_Ed[l]+E3;
}
if((Phantom_xpre[l]>=50) && (Phantom_xpre[l]<=100))
{
E4=Phantom_Ed[l]+E4;
}

}//end of steps loop

//printing out to screen the energy in each section

cout << "Energy deposited in first section = " << E1 << " MeV" << endl;
cout << "Energy deposited in second section = " << E11 << "MeV" << endl;
cout << "Energy deposited in third section = " << E12 << "MeV" << endl;
cout << "Energy deposited in fourth section = " << E13 << "MeV" << endl;
cout << "Energy deposited in fifth section = " << E14 << "MeV" << endl;
cout << "Energy deposited in second quarter = " << E2 << "MeV" << endl;
cout << "Energy deposited in third quarter = " << E3 << "MeV" << endl;
cout << "Energy deposited in fourth quarter = " << E4 << "MeV" << endl;


}//end of primary loop


NB i will be changing so that the x intervals are the same length

The Attempt at a Solution



i would assume that the best way to deal with the x-axis is, first of all define a graph (TGraph or TGraph2D which, i don't know the documentation confused me)

in my head i would write some code that first of all sets the x position to be in the middle of the interval. make an array of the E's and somehow put it into graph form.

I realize that my "attempt at a solution" isn't fantastic but I've only been coding for a couple of months and writing the above code was rather testing and now i really don't know how to proceed further. Hence any help will be greatly appreciated

Thanks

NB: I am using CERN root package
 
Physics news on Phys.org
lavster said:

Homework Statement



Im having difficulty drawing a specific graph. I am not entirely sure how to explain my question/ the graph i need to draw so i think the best way is to provide a bit of my code and explain using it the graph using it...

the code below gives a read out of the energy deposited (Phantom_Ed) for a particular x interval (Phantom_xpre). These read-outs are given a label - E1 for between -100 and -90mm, E2 for -90 and -80mm etc. What id like to do is produce a graph of the position (Phantom_xpre) with the corresponding value of Energy deposited (E1 for first interval etc) plotted on the y axis. However I am not sure how I am to go about doing this...

Homework Equations



my code at the moment (with the bit before defining the tree etc missing) is:
Your code would be easier to write if it was indented and placed inside [ code] and [ /code] tags (without the leading spaces.
lavster said:
Code: 
Long64_t nentries = T->GetEntries();// so can determine how many entries for loop
float E1=0; // defining variables for loop and setting to zero
float E2=0;
float E3=0;
float E4=0;
float E11=0;
float E12=0;
float E13=0;
float E14=0;

for(Long64_t k=0;k<nentries;k++)//original loop - gets the entry
{

  T->GetEntry(k);

  for( Long64_t l=0; l < Phantom_Nhits; l++ ) //steps loop
  {
    if((Phantom_xpre[l]>=-100) && (Phantom_xpre[l]<-90)) //if statement defining the region of x-axis of interest
    {
      E1=Phantom_Ed[l]+E1; //adding up the total energy depsoited in this region
    }

    if((Phantom_xpre[l]>=-90) && (Phantom_xpre[l]<=-80))
    {
      E11=Phantom_Ed[l]+E11;
    }

    if((Phantom_xpre[l]>=-80) && (Phantom_xpre[l]<=-70))
    {
      E12=Phantom_Ed[l]+E12;
    }

    if((Phantom_xpre[l]>=-70) && (Phantom_xpre[l]<=-60))
    {
      E13=Phantom_Ed[l]+E13;
    }

    if((Phantom_xpre[l]>=-60) && (Phantom_xpre[l]<=-50))
    {
      E14=Phantom_Ed[l]+E14;
    }

    if((Phantom_xpre[l]>=-50) && (Phantom_xpre[l]<0))
    {
      E2=Phantom_Ed[l]+E2;
    }

    if((Phantom_xpre[l]>=0) && (Phantom_xpre[l]<50))
    {
      E3=Phantom_Ed[l]+E3;
    }

    if((Phantom_xpre[l]>=50) && (Phantom_xpre[l]<=100))
    {
      E4=Phantom_Ed[l]+E4;
    }

  }//end of steps loop

  //printing out to screen the energy in each section

  cout << "Energy deposited in first section = " << E1 << " MeV" << endl;
  cout << "Energy deposited in second section = " << E11 << "MeV" << endl;
  cout << "Energy deposited in third section = " << E12 << "MeV" << endl;
  cout << "Energy deposited in fourth section = " << E13 << "MeV" << endl;
  cout << "Energy deposited in fifth section = " << E14 << "MeV" << endl;
  cout << "Energy deposited in second quarter = " << E2 << "MeV" << endl;
  cout << "Energy deposited in third quarter = " << E3 << "MeV" << endl;
  cout << "Energy deposited in fourth quarter = " << E4 << "MeV" << endl;


}//end of primary loop
This is actually a fairly simple program. One thing I noticed is that your intervals are not consistent, which could lead to some values being counted twice. For example, you have this near the beginning:
if((Phantom_xpre[l]>=-100) && (Phantom_xpre[l]<-90))

In the next if statement you have this:
if((Phantom_xpre[l]>=-90) && (Phantom_xpre[l]<=-80))
and in the following one you have this:
if((Phantom_xpre[l]>=-80) && (Phantom_xpre[l]<=-70))

If Phantom_xpre[l] is -90, the first if expression will be false and the second one will be true. If Phantom_xpre[l] is -80, both the second and third if expressions will be true, so you will be incrementing two different buckets for one value, which you don't want.

To graph this stuff you'll need to take another look at the documentation for the graphing utilities that are available to you.
lavster said:
NB i will be changing so that the x intervals are the same length

The Attempt at a Solution



i would assume that the best way to deal with the x-axis is, first of all define a graph (TGraph or TGraph2D which, i don't know the documentation confused me)

in my head i would write some code that first of all sets the x position to be in the middle of the interval. make an array of the E's and somehow put it into graph form.

I realize that my "attempt at a solution" isn't fantastic but I've only been coding for a couple of months and writing the above code was rather testing and now i really don't know how to proceed further. Hence any help will be greatly appreciated

Thanks

NB: I am using CERN root package
 
ah yes thanks :)

I think one of my main problems is ascertaining what type of graph i need - am i right in saying its a TGraph? and not TGraph2D?

thanks
 
I'm not familiar with TGraph or TGraph2D. You'll need to look at the documentation for these functions to see which is appropriate for what you're trying to do.
 
This is the code that i have now produced:

Code: 
  //defining hist

  TH1F* h = new TH1F("h", "x_{pre}", 39, -100, 100);//defining a histogram
 

  // defining variables for loop and setting to zero

  float E1=0; 
  float E2=0;
  float E3=0;
  float E4=0;
  float E5=0;
  float E6=0;
  float E7=0;
  float E8=0;
  float E9=0;
  float E10=0;
  float E11=0;
  float E12=0;
  float E13=0;
  float E14=0;
  float E15=0;
  float E16=0;
  float E17=0;
  float E18=0;
  float E19=0;
  float E20=0;
  float E21=0;
  float E22=0;
  float E23=0;
  float E24=0;
  float E25=0;
  float E26=0;  
  float E27=0;
  float E28=0;
  float E29=0;
  float E30=0;
  float E31=0;
  float E32=0;
  float E33=0;
  float E34=0;
  float E35=0;
  float E36=0;
  float E37=0;
  float E38=0;
  float E39=0;
  
  
  //loop to find energy deposited in each x interval

  Long64_t nentries = T->GetEntries();// so can determine how many entries for loop
  
  for(Long64_t k=0;k<nentries;k++)//original loop - gets the entry
    {
      T->GetEntry(k);
      for( Long64_t l=0; l < Phantom_Nhits; l++ ) //steps loop
	{
	  if((Phantom_xpre[l]>=-100) && (Phantom_xpre[l]<-95)) //if statement defining the region of x-axis of interest
	    {
	      E1=Phantom_Ed[l]+E1; //adding up the total energy depsoited in this region
	    }
	  
	  if((Phantom_xpre[l]>=-95) && (Phantom_xpre[l]<-90))
	    {
	      E2=Phantom_Ed[l]+E2;
	    }
	  if((Phantom_xpre[l]>=-90) && (Phantom_xpre[l]<-85))
	    {
	      E3=Phantom_Ed[l]+E3;
	    }
	  if((Phantom_xpre[l]>=-85) && (Phantom_xpre[l]<-80))
	    {
	      E4=Phantom_Ed[l]+E4;
	    }
	  
	  if((Phantom_xpre[l]>=-80) && (Phantom_xpre[l]<-75))
	    {
	      E5=Phantom_Ed[l]+E5;
	    }
	  
	  if((Phantom_xpre[l]>=-75) && (Phantom_xpre[l]<-70))
	    {
	      E6=Phantom_Ed[l]+E6;
	    }
	  
	  if((Phantom_xpre[l]>=-75) && (Phantom_xpre[l]<-65))
	    {
	      E7=Phantom_Ed[l]+E7;
	    }
	  
	  if((Phantom_xpre[l]>=-65) && (Phantom_xpre[l]<-60))
	    {
	      E8=Phantom_Ed[l]+E8;
	    }
	  
	  if((Phantom_xpre[l]>=-60) && (Phantom_xpre[l]<-55))
	    {
	      E9=Phantom_Ed[l]+E9;
	    }
	  
	  if((Phantom_xpre[l]>=-55) && (Phantom_xpre[l]<-50))
	    {
	      E10=Phantom_Ed[l]+E10;
	    }
	  
	  if((Phantom_xpre[l]>=-50) && (Phantom_xpre[l]<-45))
	    {
	      E11=Phantom_Ed[l]+E11;
	    }
	  
	  if((Phantom_xpre[l]>=-45) && (Phantom_xpre[l]<-40))
	    {
	      E12=Phantom_Ed[l]+E12;
	    }
	  
	  if((Phantom_xpre[l]>=-40) && (Phantom_xpre[l]<-35))
	    {
	      E13=Phantom_Ed[l]+E13;
	    }
	  
	  
	  if((Phantom_xpre[l]>=-35) && (Phantom_xpre[l]<-30))
	    {
	      E14=Phantom_Ed[l]+E14;
	    }
	  
	  if((Phantom_xpre[l]>=-30) && (Phantom_xpre[l]<-25))
	    {
	      E15=Phantom_Ed[l]+E15;
	    }
	  
	  if((Phantom_xpre[l]>=-25) && (Phantom_xpre[l]<-20))
	    {
	      E16=Phantom_Ed[l]+E16;
	    }
	  
	  if((Phantom_xpre[l]>=-20) && (Phantom_xpre[l]<-15))
	    {
	      E17=Phantom_Ed[l]+E17;
	    }
	  
	  if((Phantom_xpre[l]>=-15) && (Phantom_xpre[l]<-10))
	    {
	      E18=Phantom_Ed[l]+E18;
	    }
	  
	  
	  if((Phantom_xpre[l]>=-10) && (Phantom_xpre[l]<-5))
	    {
	      E19=Phantom_Ed[l]+E19;
	    }
	  
	  if((Phantom_xpre[l]>=-5) && (Phantom_xpre[l]<0))
	    {
	      E20=Phantom_Ed[l]+E20;
	    }
	  
	  if((Phantom_xpre[l]>=0) && (Phantom_xpre[l]<5))
	    {
	      E21=Phantom_Ed[l]+E21;
	    }
	  
	  if((Phantom_xpre[l]>=5) && (Phantom_xpre[l]<10))
	    {
	      E22=Phantom_Ed[l]+E22;

	    }
	  
	  if((Phantom_xpre[l]>=10) && (Phantom_xpre[l]<15))
	    {
	      E23=Phantom_Ed[l]+E23;
	    }
	  
	  
	  if((Phantom_xpre[l]>=15) && (Phantom_xpre[l]<20))
	    {
	      E23=Phantom_Ed[l]+E23;
	    }
	  if((Phantom_xpre[l]>=20) && (Phantom_xpre[l]<25))
	    {
	      E24=Phantom_Ed[l]+E24;
	    }
	  
	  if((Phantom_xpre[l]>=25) && (Phantom_xpre[l]<30))
	    {
	      E25=Phantom_Ed[l]+E25;
	    }
	  if((Phantom_xpre[l]>=30) && (Phantom_xpre[l]<35))
	    {
	      E26=Phantom_Ed[l]+E26;
	    }
	  if((Phantom_xpre[l]>=35) && (Phantom_xpre[l]<40))
	    {
	      E27=Phantom_Ed[l]+E27;
	    }
	  if((Phantom_xpre[l]>=40) && (Phantom_xpre[l]<45))
	    {
	      E28=Phantom_Ed[l]+E28;
	    }
	  if((Phantom_xpre[l]>=45) && (Phantom_xpre[l]<50))
	    {
	      E29=Phantom_Ed[l]+E29;
	    }
	  if((Phantom_xpre[l]>=50) && (Phantom_xpre[l]<55))
	    {
	      E30=Phantom_Ed[l]+E30;
	    }
	  if((Phantom_xpre[l]>=55) && (Phantom_xpre[l]<60))
	    {
	      E31=Phantom_Ed[l]+E31;
	    }
	  if((Phantom_xpre[l]>=60) && (Phantom_xpre[l]<65))
	    {
	      E32=Phantom_Ed[l]+E32;
	    }
	  if((Phantom_xpre[l]>=65) && (Phantom_xpre[l]<70))
	    {
	      E33=Phantom_Ed[l]+E33;
	    }
	  if((Phantom_xpre[l]>=70) && (Phantom_xpre[l]<75))
	    {
	      E34=Phantom_Ed[l]+E34;
	    }
	  if((Phantom_xpre[l]>=75) && (Phantom_xpre[l]<80))
	    {
	      E35=Phantom_Ed[l]+E35;
	    }
	  if((Phantom_xpre[l]>=80) && (Phantom_xpre[l]<85))
	    {
	      E36=Phantom_Ed[l]+E36;
	    }
	  if((Phantom_xpre[l]>=85) && (Phantom_xpre[l]<90))
	    {
	      E37=Phantom_Ed[l]+E37;
	    }
	  if((Phantom_xpre[l]>=90) && (Phantom_xpre[l]<95))
	    {
	      E38=Phantom_Ed[l]+E38;
	    }
	  if((Phantom_xpre[l]>=95) && (Phantom_xpre[l]<=100))
	    {
	      E39=Phantom_Ed[l]+E39;
	    }
	  h->Fill(Phantom_xpre[k]);
	}//end of steps loop
      
    }//end of primary loop
  //printing out to screen the energy in each section
  
  cout << "Energy deposited in 1 section = " << E1 << "\tMeV" << endl;
  cout << "Energy deposited in 2 section = " << E2 << "MeV" << endl;
  cout << "Energy deposited in 3 section = " << E3 << "MeV" << endl;
  cout << "Energy deposited in 4 section = " << E4 << "MeV" << endl;
  cout << "Energy deposited in 5 section = " << E5 << "MeV" << endl;
  cout << "Energy deposited in 6 section = " << E6 << "MeV" << endl;
  cout << "Energy deposited in 7 section = " << E7 << "MeV" << endl;
  cout << "Energy deposited in 8 section = " << E8 << "MeV" << endl;
  cout << "Energy deposited in 9 section = " << E9 << "MeV" << endl;
  cout << "Energy deposited in 10 section = " << E10 << "MeV" << endl;
  cout << "Energy deposited in 11 section = " << E11 << "MeV" << endl;
  cout << "Energy deposited in 12 section = " << E12 << "MeV" << endl;
  cout << "Energy deposited in 13 section = " << E13 << "MeV" << endl;
  cout << "Energy deposited in 14 section = " << E14 << "MeV" << endl;
  cout << "Energy deposited in 15 section = " << E15 << "MeV" << endl;
  cout << "Energy deposited in 16 section = " << E16 << "MeV" << endl;
  cout << "Energy deposited in 17 section = " << E17 << "MeV" << endl;
  cout << "Energy deposited in 18 section = " << E18 << "MeV" << endl;
  cout << "Energy deposited in 19 section = " << E19 << "MeV" << endl;
  cout << "Energy deposited in 20 section = " << E20 << "MeV" << endl;
  cout << "Energy deposited in 21 section = " << E21 << "MeV" << endl;
  cout << "Energy deposited in 22 section = " << E22 << "MeV" << endl;
  cout << "Energy deposited in 23 section = " << E23 << "MeV" << endl;
  cout << "Energy deposited in 24 section = " << E24 << "MeV" << endl;
  cout << "Energy deposited in 25 section = " << E25 << "MeV" << endl;
  cout << "Energy deposited in 26 section = " << E26 << "MeV" << endl;
  cout << "Energy deposited in 27 section = " << E27 << "MeV" << endl;
  cout << "Energy deposited in 28 section = " << E28 << "MeV" << endl;
  cout << "Energy deposited in 29 section = " << E29 << "MeV" << endl;
  cout << "Energy deposited in 30 section = " << E30 << "MeV" << endl;
  cout << "Energy deposited in 31 section = " << E31 << "MeV" << endl;
  cout << "Energy deposited in 32 section = " << E32 << "MeV" << endl;
  cout << "Energy deposited in 33 section = " << E33 << "MeV" << endl;
  cout << "Energy deposited in 34 section = " << E34 << "MeV" << endl;
  cout << "Energy deposited in 35 section = " << E35 << "MeV" << endl;
  cout << "Energy deposited in 36 section = " << E36 << "MeV" << endl;
  cout << "Energy deposited in 37 section = " << E37 << "MeV" << endl;
  cout << "Energy deposited in 38 section = " << E38 << "MeV" << endl;
  cout << "Energy deposited in 39 section = " << E39 << "MeV" << endl;

  //setting up an array

  Float_t Es []={0, E1, E2, E3, E4, E5, E6, E7, E8, E9, E10, E11, E12, E13, E14, E15, E16, E17, E18, E19, E20, E21, E22, E23, E24, E25, E26, E27, E28, E29, E30, E31, E32, E33, E34, E35, E36, E37, E38, E39};// array with 0 at beginning because of for loop below-bin zero is an underflow bin so m>=1.

  //loop for the graph
  
cout << "debug " << Es[2] << "MeV" << endl;

  
  // Float_t x[10000];
  //Float_t y[10000]; //defining the x and y variables (arrays)
  
 

cout << "bonjour" << endl;
 
 TGraph *g =new TGraph(h->GetNbinsX());//defining a graph
 // g->SetTitle("hi");

 cout << "hello" << endl;
 float graphPointNumber=0;//where we are
  cout << "salut" << endl;
  for(int m=1;  m<=h->GetNbinsX();m++)
   {
      cout << "hi" << endl;
      cout << "m="<< m  << endl;

     g->SetPoint(graphPointNumber,h->GetBinCenter(m),Es[m]); 
     graphPointNumber++;
 cout << "value =" << Es[m] << "MeV" << endl;

   }
 
 //drawing the graph
 TCanvas* c = new TCanvas("c","I am a graph",1000,1000);
 c->Divide(1,1);
 c->cd(1);
 g->SetTitle("hilhhhhkggggggggggggggggggggggggggggg");
 g->SetLineWidth(4);
 g->SetLineColor(4);
 g->GetXaxis()->SetTitle("x position (mm)");
 g->GetYaxis()->SetTitle("Energy Deposited (MeV)");
 // g->SetName("hope this works");
 g->Draw("ACP"); 
 c->cd();

however, it is really long and messy. is there a way i could make it neater? like have loops over the x direction perhaps? and also, my title is not appearing and i don't no what I am doing wrong...

thanks
 
I don't see anything wrong with the long list of if statements, but you could get rid of the 39 En variables by storing the values in an array with 39 elements.

I don't know what you mean when you say "loops over the x direction."

Are any titles appearing? It looks like you have titles for the x-axis and y-axis. Are these appearing? What about this title - "hilhhhhkggggggggggggggggggggggggggggg"?

What do the statements c->cd(1) and c->cd() do?
 
ah kl thanks :)

yes, the axis titles are working but not the overall graph title, which is what that's meant to be.

and c->cd(1) and c->cd() arent really needed i don't think.the first one says go to pallete one on canvas, but there is only one, and the c->cd() was just me trying to see if it makes a difference. i didnt have it in before...

thanks
 
Since what you are using seems to be a specialized graphics package, you need to look at the existing documentation for that package, particularly the SetTitle function on the TGraph class/struct.

One other thing I noticed is in this code:
Code: 
for(int m=1;  m<=h->GetNbinsX();m++)
{
   cout << "hi" << endl;
   cout << "m="<< m  << endl;

   g->SetPoint(graphPointNumber,h->GetBinCenter(m),Es[m]); 
   graphPointNumber++;
   cout << "value =" << Es[m] << "MeV" << endl;
}
In this loop you are checking whether m <= h->GetNbinsX() each iteration of the loop. If GetNbinsX() returns the same value every time, you could do this:
Code: 
int nBins = h -> GetNbinsX();
for(int m=1;  m <= nBins; m++)
{
   cout << "hi" << endl;
   cout << "m="<< m  << endl;

   g->SetPoint(graphPointNumber,h->GetBinCenter(m),Es[m]); 
   graphPointNumber++;
   cout << "value =" << Es[m] << "MeV" << endl;
}
I would also get rid of the cout statements once this is working correctly.

Also, I don't know what the type of the first argument to SetPoint should be, but I suspect it should be int or long, and you have declared graphPointNumber to be of type float. If you are counting things you shouldn't be using float or double.
 
Mark44 said:
I don't see anything wrong with the long list of if statements, but you could get rid of the 39 En variables by storing the values in an array with 39 elements.
I don't know what you mean when you say "loops over the x direction."
i mean having a loop so that instead of writing it all out we just right one if statement for x_pre by eg defining a variable x=0 and so the interval would be between x and x+10 and then the next interation would be between this x+10 (say equal to y) and y+10 and so on until the end. And instead of specifing Es1 etc it would just be Es. but I am not sure how to do that, or if it is at all possible...

this is the kind of thing i was thinking but it won't work but i can't get my head around it.

Code: 
int x=-100;

for(long64_t=0; l<Phantom_Nhits; l++)
{
if((Phantom_xpre[l]>=x)&&(Phantom_xpre[l]<x+10)
{
E[l]=Phantom_Ed[l]+E[l]
}
}

but this won't do what i wanted it to do (ie its not consistent with what i got before...
im goin to be expanding to three dimensions and it will get very confusing otherwise...
and by this array do you mean have Es[]={Es1, Es2,...Es39}={0} at the beginning. how about defining it to be a float or whatever?
Mark44 said:
Also, I don't know what the type of the first argument to SetPoint should be, but I suspect it should be int or long, and you have declared graphPointNumber to be of type float. If you are counting things you shouldn't be using float or double.
g is a graph. so for counting things it should be an int? also using eg graphPointNumber++ increases it by one...is there anyway to increase by eg ten by eg putting more ++ at the end?
 
  • #10
lavster said:
i mean having a loop so that instead of writing it all out we just right one if statement for x_pre by eg defining a variable x=0 and so the interval would be between x and x+10 and then the next interation would be between this x+10 (say equal to y) and y+10 and so on until the end. And instead of specifing Es1 etc it would just be Es. but I am not sure how to do that, or if it is at all possible...

this is the kind of thing i was thinking but it won't work but i can't get my head around it.
You could do this with a for loop nested inside another for loop. The inner for loop would work with the values within an interval [x, x + 10]. The outer loop would take care of the four intervals.

Code: 
for (int i = 0; i < 40; i += 10)
{
   for (int j = 0; j < 10; j++)
   {
      // Do stuff
   }
}
I'm going to indent your code below to make it more readable.
lavster said:
Code: 
int x=-100;

for(long64_t=0; l<Phantom_Nhits; l++)
{
   if((Phantom_xpre[l]>=x)&&(Phantom_xpre[l]<x+10)
   {
      E[l]=Phantom_Ed[l]+E[l]
   }
}
Comments on the code above.
1. The first statement in the loop header is wrong. long64_t is apparently a type. The expression should long64_t l = 0;
2. I don't see any need for a 64-bit loop counter variable when the array has only 40 elements, unless the default size for an int in your machine architecture is 64 bits. An int will probably do just fine.
lavster said:
but this won't do what i wanted it to do (ie its not consistent with what i got before...
im goin to be expanding to three dimensions and it will get very confusing otherwise...
and by this array do you mean have Es[]={Es1, Es2,...Es39}={0} at the beginning. how about defining it to be a float or whatever?
Float or double, your call.

The code below will declare and define a 40-element array of type double, and then set all 40 elements to 0.0. Note that the array elements run from Es[0] through Es[39]. Arrays in C/C++ are zero-based, meaning the first element has an index of 0.
Code: 
double Es[40];

for (i = 0; i < 40; i++)
{
   Es[i] = 0.0;
}
lavster said:
g is a graph. so for counting things it should be an int? also using eg graphPointNumber++ increases it by one...is there anyway to increase by eg ten by eg putting more ++ at the end?
If you're counting things, yes, you should use an int or some other integral type, depending on how many things you need to count.

The ++ operator increments a variable by 1. There is no operator for incrementing by larger amounts, but there is an operator that adds whatever to a variable and stores the new value in the variable, +=.

Code: 
count += 10;
This is functionally equivalent to writing count = count + 10;
 
  • #11
Thank you very much for all your help! However i have another question...

in the following bit of code:


Code: 
for (int i = -100; i < 100; i += 5)
{

\\ code which means "if(the value of Phantom_xpre[] lies between the value of i and the next increment (ie here i and i+5) then" (so eg if phantom_xpre was -100 then it would be in the 1st interval, if it was -90 it would be in the second and so on...)

   for (int j = 0; j < 40; j++)
   {
      Es[i+100]=Es[i+100]+Phantom_Ed[j]
   }
}

This is meant to say if the value of Phantom_xpre is between i and i+5 then add the value of Phantom_Ed to the original value of Es for this interval. so if my kth entry , phantom_xpre=-90 then it will be in the third interval so the corresponding value of Phantom_Ed will be added to the existing total of of the energy deposited (Phantom_Ed) in interval 3 (ie it will be added to E[3]. does this make sense? is this what my code is doing? will i get 40 different values of Es? and what should be the code that replaces the bit in the inverted commas? is this basically doing what i did before with my previous code?

thanks
 
  • #12
this is what i have now. it does not work :(

Code: 
  // defining variables for loop and setting to zero
   
   double Es[40];
   
   for (int i=0; i<40;i++)
     {
       Es[i]=0.0;
     }
   
   
   //loop to find energy deposited in each x interval
   
   Long64_t nentries = T->GetEntries();// so can determine how many entries for loop
   
   for(Long64_t k=0;k<nentries;k++)//original loop - gets the entry
     {
       T->GetEntry(k);
       for( Long64_t l=0; l < Phantom_Nhits; l++ ) //steps loop
	 {
	   for(int m=0;m<195;m+=5)
	     {
	       if((Phantom_xpre[l]>=(m-100)) && (Phantom_xpre[m]<(m-95))) //if statement defining the region of x-axis of interest
		 {
		   // for(int n=0;n<5;n++)
		     //{
		       Es[m/5]=Phantom_Ed[l]+Es[m/5]; //adding up the total energy depsoited in this region
		       //  }
		 }
	     }
	   
	   h->Fill(Phantom_xpre[k]);
	 }//end of steps loop
       
     }//end of primary loop
   
   //printing out to screen the energy in each section
  
  for (int i=0; i<40;i++)
     {
       cout << "Energy deposited in" << i << "  section = " << Es[i] << "MeV" << endl;
     }
   
 cout << "bonjour" << endl;
 
 TGraph *g =new TGraph(h->GetNbinsX());//defining a graph
 g->SetTitle("hi");

 // cout << "hello" << endl;
 int graphPointNumber=0;//where we are
 // cout << "salut" << endl;
  for(int p=0;  p<=h->GetNbinsX();p++)
   {
      cout << "hi" << endl;
      cout << "m="<< p  << endl;

     g->SetPoint(graphPointNumber,h->GetBinCenter(p+1),Es[p]); 
     graphPointNumber++;
     cout << "value =" << Es[p] << "MeV" << endl;

   }
 
 //drawing the graph
 TCanvas* c = new TCanvas("c","I am a graph",1000,1000);
 g->SetTitle("electrons 25 MeV");
 g->SetLineWidth(4);
 g->SetLineColor(4);
 g->GetXaxis()->SetTitle("x position (mm)");
 g->GetYaxis()->SetTitle("Energy Deposited (MeV)");
 g->Draw("ACP"); 
 c->cd();
 
  • #13
lavster said:
Thank you very much for all your help! However i have another question...

in the following bit of code:
Your comment was so long that it extended way beyond the display window. I have split it into several shorter lines. Also, the C++ comment delimiter is //, not \\.
lavster said:
Code: 
for (int i = -100; i < 100; i += 5)
{

// code which means "if(the value of Phantom_xpre[] lies between the value
// of i and the next increment (ie here i and i+5) then" (so eg if phantom_xpre
// was -100 then it would be in the 1st interval, if it was -90 it would be in
// the second and so on...)

   for (int j = 0; j < 40; j++)
   {
      Es[i+100]=Es[i+100]+Phantom_Ed[j]
   }
}

This is meant to say if the value of Phantom_xpre is between i and i+5 then add the value of Phantom_Ed to the original value of Es for this interval. so if my kth entry , phantom_xpre=-90 then it will be in the third interval so the corresponding value of Phantom_Ed will be added to the existing total of of the energy deposited (Phantom_Ed) in interval 3 (ie it will be added to E[3]. does this make sense? is this what my code is doing? will i get 40 different values of Es? and what should be the code that replaces the bit in the inverted commas? is this basically doing what i did before with my previous code?

thanks

This won't work correctly. Presumably you have declared Es as an array with 40 elements, so you definitely don't want to be putting things in Es[i + 100]. Think about it: when i is 95, you are going to be trying to assign a value to Es[195], which is waaaaaaay beyond of the array you declared, and will cause bad (or even terrible) things to happen.

There are a couple of ways to approach this, but I'll discuss only one. Your loop counter runs by 5s from -100 to 95, inclusive, with increments of 5.

You want the loop indexes to be paired with array indexes in this manner:
Code: 
Loop   Array index
-100   0
-95     1
-90     2
-85     3
.
.
.
-5     19
0     20
5     21
10     22
.
.
.
95     39
To get this pairing, add 100 to the array index and divide by 5. For example, if the loop index is -85, the array index is (-85 + 100)/5 = 15/5 = 3.

If the loop index is 95, the array index is (95 + 100)/5 = 195/5 = 39.
 
  • #14
Note: I modified the code below by taking out the lines of code you commented. I also tightened up the indenting so that the lines of code didn't extend so far to the right.
lavster said:
Code: 
Long64_t nentries = T->GetEntries();// so can determine how many entries for loop

for(Long64_t k=0;k<nentries;k++)//original loop - gets the entry
{
   T->GetEntry(k);
   for( Long64_t l=0; l < Phantom_Nhits; l++ ) //steps loop
   {
      for(int m=0;m<195;m+=5)
      {
         if((Phantom_xpre[l]>=(m-100)) && (Phantom_xpre[m]<(m-95))) //if statement defining the region of x-axis of interest
         {
            Es[m/5]=Phantom_Ed[l]+Es[m/5]; //adding up the total energy depsoited in this region
         }
      }

      h->Fill(Phantom_xpre[k]);
   }//end of steps loop
       
}//end of primary loop
First off, your loops are nested too deeply at three deep. If I understand what you're trying to do, it can be done one loop nested inside another.

This is how I understand things. The first array has lots of raw data, maybe thousands or millions (or whatever) of meV readings. What you need to do is to put each of them into one of 40 buckets so that you can get a histogram of this data.

Apparently most of the data is between -100 (somethings) and +100 (somethings). Let's call the raw data x_pre for short.
If -100 <= x_pre && x_pre < -95, then increment bucket[0].
If -95 <= x_pre && x_pre < -90, then increment bucket[1].
And so on.

The code would look like this:
Code: 
for (i = 0; i < nData; i++)
{
   for(j = -100; j < 100; j += 5)
   {
      if ( j <= x_pre[i] && x_pre[i] < j + 5)
      {
         bucket[(j + 100)/5]++;
      }
   }
}

Out of curiosity, what does T -> GetEntry(k); do? If it's supposed to return a value, you aren't storing the value anywhere, so this is just wasted effort in this case.
Also, what does h->Fill(Phantom_xpre[k]); do?
 
  • #15
okay... i think i am getting it slowly but surely...

This is whast I am trying to do in words:

I am modelling a beam path of eg protons in a phantom (cylinder of water) each particle that gets fired into the cylinder is an event k, and the MC program follows each k event and simulates how it interacts with the phantom. each time it interacts, step l, it produces more particles. The primary proton k and each particle produced in the step l can deposit energy. I want to produce a graph which gives a measure of how the energy deposted varies with depth. In order to do this, i need to split up the x distance (Phantom_xpre) into 5mm intervals and then record the energy deposted (Phantom_Ed) in each of these intervals.

so T->GetEntry[k] gets the event that I am interested in and follows its path.
for me to have a graph with the x-axis being Phantom_xpre and y-axis being the corresponding total energy deposited in the interval i thought the easiest way was to define a histogram for Phantom_xpre and then fill it (h->Fill(Phantom_xpre)) and then when drawng the graph the x value that corresponds to the energy deposited was the middle of the bin in the histogram

so yes you have the idea. but it looks to me that its just counting the number of entries that are in each interval? so instead of "bucket[(j + 100)/5]++;" i need code that adds the Phantom_Ed to the total Phantom_Ed measured in each bucket? i thought this would that be something like: bucket[(j+100)/5]=Phantom_Ed+bucket[(j+100/5)]? or am i getting the definition of "bucket[(j + 100)/5]++;" wrong? also what I am saying is wrong as it gets really bizarre values

thanks
 
  • #16
lavster said:
okay... i think i am getting it slowly but surely...

This is whast I am trying to do in words:

I am modelling a beam path of eg protons in a phantom (cylinder of water)
This might be a language translation problem. In English a phantom is something apparently seen, heard, or sensed, but having no physical reality; a ghost or an apparition. A cylinder of water would be a column of water.
lavster said:
each particle that gets fired into the cylinder is an event k, and the MC program follows each k event and simulates how it interacts with the phantom. each time it interacts, step l, it produces more particles. The primary proton k and each particle produced in the step l can deposit energy. I want to produce a graph which gives a measure of how the energy deposted varies with depth. In order to do this, i need to split up the x distance (Phantom_xpre) into 5mm intervals and then record the energy deposted (Phantom_Ed) in each of these intervals.
Your explanation is somewhat confusing. When you say "an event k" or a "primary proton k" that suggests that events and protons are numbers. It makes more sense to talk about the k-th event or k-th primary proton.
lavster said:
so T->GetEntry[k] gets the event that I am interested in and follows its path.
In your code you have T -> GetEntry(k), which suggests that GetEntry() is a function, and probably one that returns a value of some kind. If that's the case, you are using this function incorrectly, since you are not storing the return value anywhere.

One way to store the return value is to assign the value returned to a variable, like this:
result = T -> GetEntry(k);

lavster said:
for me to have a graph with the x-axis being Phantom_xpre and y-axis being the corresponding total energy deposited in the interval i thought the easiest way was to define a histogram for Phantom_xpre and then fill it (h->Fill(Phantom_xpre)) and then when drawng the graph the x value that corresponds to the energy deposited was the middle of the bin in the histogram

so yes you have the idea. but it looks to me that its just counting the number of entries that are in each interval? so instead of "bucket[(j + 100)/5]++;" i need code that adds the Phantom_Ed to the total Phantom_Ed measured in each bucket? i thought this would that be something like: bucket[(j+100)/5]=Phantom_Ed+bucket[(j+100/5)]? or am i getting the definition of "bucket[(j + 100)/5]++;" wrong? also what I am saying is wrong as it gets really bizarre values

It depends on what you're putting in the buckets. If you are just counting the number of entries in each interval, then when you find an entry in a given interval, increment the bucket for that interval. That's exactly what my code does -- bucket[(j + 100)/5]++;

On the other hand, if the buckets contain the total energy of the entries in the interval, then you should increment the bucket for that interval by the energy for that entry, which you can do in either of two ways:
bucket[(j+100)/5] = bucket[(j+100/5)] + Phantom_Ed;
(The above is similar to how you wrote it.)
OR
bucket[(j+100)/5] += Phantom_Ed;
These two ways do the same thing.

By the way, Phantom_Ed does not seem to me to be a useful variable name. Unless you can convince me that there is some significance to the "Phantom" part, I would get rid of it in this and the other variables that have this prefix. The important part of the name is E, which is energy. I guess "d" signifies distance, but I'm not sure.

Same applies to Phantom_xpre. The "x" part indicates distance (maybe depth in the water column?), but I have no idea what the "pre" part is supposed to mean. Well-chosen variable names should convey to a reader exactly what they are being used for.
 
  • #17
Mark44 said:
This might be a language translation problem. In English a phantom is something apparently seen, heard, or sensed, but having no physical reality; a ghost or an apparition. A cylinder of water would be a column of water.

In radiotherapy, which is what my project is on a phantom is "an object generally comprising of tissue substitute materials used to simulate a patient or part thereof".

Your explanation is somewhat confusing. When you say "an event k" or a "primary proton k" that suggests that events and protons are numbers. It makes more sense to talk about the k-th event or k-th primary proton.

okay thanks :) (im not very good with this terminology at all lol)

It depends on what you're putting in the buckets. If you are just counting the number of entries in each interval, then when you find an entry in a given interval, increment the bucket for that interval. That's exactly what my code does -- bucket[(j + 100)/5]++;

On the other hand, if the buckets contain the total energy of the entries in the interval, then you should increment the bucket for that interval by the energy for that entry, which you can do in either of two ways:
bucket[(j+100)/5] = bucket[(j+100/5)] + Phantom_Ed;
(The above is similar to how you wrote it.)
OR
bucket[(j+100)/5] += Phantom_Ed;
These two ways do the same thing.


ah excellent- thanks!

By the way, Phantom_Ed does not seem to me to be a useful variable name. Unless you can convince me that there is some significance to the "Phantom" part, I would get rid of it in this and the other variables that have this prefix. The important part of the name is E, which is energy. I guess "d" signifies distance, but I'm not sure.

Same applies to Phantom_xpre. The "x" part indicates distance (maybe depth in the water column?), but I have no idea what the "pre" part is supposed to mean. Well-chosen variable names should convey to a reader exactly what they are being used for.

These were the variables that were given to me - the d is deposited and this is concerned with what is happening in the phantom. there are other variables called Prim that is to do with something else and it is to distinguish between the two. yes, x is depth into cylinder, when its lying on its curvey side, pre is to distinguish between the post position. don't no what that is but I am using someone elses variables.

Thank you very much for all your help! It is much appreciated. :D
 

Similar threads

Comp Sci C++ how to obtain the x value of graph corresponding to a certain y value
Replies
5
Views
3K
Comp Sci Solve C++ & ROOT Problems: 3 Questions Answered
Replies
8
Views
2K
Comp Sci Matrix inverse with LU decomposition in C++
Replies
7
Views
8K
Comp Sci C++ (ROOT) Form a matrix and send it to a 2d Histogram
Replies
13
Views
3K
Comp Sci C++: Method of Minimum Squares used for System Solving
Replies
13
Views
2K
Comp Sci Understanding Newton's Method in C++
Replies
8
Views
3K
Comp Sci Problem with my C++ code (MATH/Physics based)
Replies
12
Views
2K
Comp Sci Calculating GPA with Loops in C++
Replies
6
Views
4K
Comp Sci C++ object oriented, classes/methods
Replies
9
Views
3K
Comp Sci Why Is My C++ Program Outputting Zeroes Instead of Averages?
Replies
5
Views
4K

Hot Threads

Recent Insights

Back
Top