Lego Plot Spherical Coordinates

From: Steve Aplin (aplin@mail.desy.de)
Date: Tue Feb 22 2000 - 16:07:50 MET


Hi 

I'm trying to display a TH2F as a lego plot using spherical coordinates.
The two dimensions of my histogram are theta and phi and with field
strength E given as the value. When I plot this in spherical coordinates
the value of r is fixed and independent of E, this gives the appearance of
a value for E where in fact the value of E is zero. Is there any option to
set the value of r to zero or a negligible value?

cheers

steve.



{

/***********************************************************************
									
Macro Name: E_field_ntuple.cxx
									
macro to create a ntuple with values taken from file: Efield.dat
  
author: Steve Aplin		Date: 
  
**********************************************************************/


gROOT->Reset("a");

Float_t c_0 = 2.99e8;		// speed of light in free space

Int_t Nv = 4;                             // Number of points of v
Int_t Nr = 10;       			  // Number of points of r
Int_t Ntheta = 9;			  // Number of points of theta
Int_t Nphi = 16;			  // Number of points of phi

Float_t v;	                          // Velocity of point charge

Float_t r;				  // radius from charge in meters
Float_t phi;				  // phi in degrees
Float_t theta;			          // theta in degrees
Float_t theta_rad;

Float_t x;                                // cartesian co-ordinate
Float_t y;                                // cartesian co-ordinate
Float_t z;                                // cartesian co-ordinate

Float_t Er; 			       	  // Value of Electric field strength 


// create File "rootfile.root" to write output to
TFile *rootfile = new TFile("rootfile.root","RECREATE");

// create Canvas "Electric Field Plots"  // ,(top right),(size hz,vrt) 
c1 = new TCanvas("c1","Electric Field Plots", 200,10,900,700);

                                  //,(bot_left),(top_right),(colour)
pad1 = new TPad("pad1","This is pad1",0.02,0.52,0.48,0.98,21);
pad2 = new TPad("pad2","This is pad2",0.52,0.52,0.98,0.98,21);
pad3 = new TPad("pad3","This is pad3",0.02,0.02,0.48,0.48,21);
pad4 = new TPad("pad4","This is pad4",0.52,0.02,0.98,0.48,21);

pad1->Draw();
pad2->Draw();
pad3->Draw();
pad4->Draw();

// create ntuple "ntuple" 
TNtuple *ntuple = new TNtuple("ntuple","E field","Er:v:r:theta:phi:x:y:z"); 

// create Histogram "Field_Map_v_0c"
TH2F *Field_Map_v_0c = new TH2F("Field_Map_v_0c",
			    "v = 0",16,0,360,10,0,1 );

// create Histogram "Field_Map_v_05c"
TH2F *Field_Map_v_05c = new TH2F("Field_Map_v_05c",
			    "v = 0.5c",16,0,360,10,0,1 );

// create Histogram "Field_Map_v_09c"
TH2F *Field_Map_v_09c = new TH2F("Field_Map_v_09c",
			    "v = 0.9c",16,0,360,10,0,1 );

// create Histogram "Field_Map_v_099c"
TH2F *Field_Map_v_099c = new TH2F("Field_Map_v_099c",
			    "v = 0.99c",16,0,360,10,0,1 );


ifstream data_file("E_field.dat");	// The input file

if (data_file.bad())			// check for bad file open
  
{    
  cerr << "Error: could not open E_field.dat\n";
  exit (8);  
}

for (Int_t i = 0; i < Nv; i++)
{
  for (Int_t j = 0; j < Nr; j++)
    {
      for (Int_t k = 0; k < Ntheta; k++)
	{
	  for (Int_t l = 0; l < Nphi; l++)
	    {
	      
	      data_file >> Er >> v >> r >> theta >> phi >> x >> y >> z;
	      //Take value from  
	      // file and assign  
	      
	      theta_rad = theta / 57.3;	      
	      
	      ntuple->Fill(Er,v,r,theta,phi,x,y,z);
	      // fill ntuple with 
	      // values 	
	      
	      if ((theta == 0) || (theta == 180)) l = Nphi;   
	      // to prevent looping over phi at singular points of theta

	    }

	  // fill histogram bin 
	  if ( i == 0 ){	 
	  Field_Map_v_0c->SetBinContent
	    (Field_Map_v_0c->GetBin((theta/180*8)+0.5, r*10-.5), Er);
	  
	  Field_Map_v_0c->SetBinContent
	    (Field_Map_v_0c->GetBin((theta/180*8)+8.5, r*10-.5), Er);

	  }
 	  
	  // fill histogram bin 
	  if ( i == 1 ){
	    Field_Map_v_05c->SetBinContent
	      (Field_Map_v_05c->GetBin((theta/180*8)+.5,(r*10)-.5), Er);
	    
	    Field_Map_v_05c->SetBinContent
	      (Field_Map_v_05c->GetBin((theta/180*8)+8.5,(r*10)-.5), Er);

	  }
	  
	  // fill histogram bin 
	  if ( i == 2 ){	 
	    Field_Map_v_09c->SetBinContent
	      (Field_Map_v_09c->GetBin((theta/180*8)+.5,(r*10)-.5), Er);
	   
	    Field_Map_v_09c->SetBinContent
	      (Field_Map_v_09c->GetBin((theta/180*8)+8.5,(r*10)-.5), Er);
	  }
	  
	  // fill histogram bin 
	  if ( i == 3 ){	 
	    Field_Map_v_099c->SetBinContent
	      (Field_Map_v_099c->GetBin((theta/180*8)+.5,(r*10)-.5), Er);

	    Field_Map_v_099c->SetBinContent
	      (Field_Map_v_099c->GetBin((theta/180*8)+8.5,(r*10)-.5), Er);  
}

	}
    }
} 

pad1->cd();
Field_Map_v_0c->Draw("lego1 pol");  

pad2->cd();
Field_Map_v_05c->Draw("lego1 pol");  

pad3->cd();
Field_Map_v_09c->Draw("lego1 pol");  

pad4->cd();
Field_Map_v_099c->Draw("lego1 pol");  


ntuple->Scan("Er:v", "(r == 1.0)&&(theta == 90)&&(phi==0)");


rootfile->Write();				// write ntuple to
						// file 

//rootfile->Close();		  


}



This archive was generated by hypermail 2b29 : Tue Jan 02 2001 - 11:50:19 MET