Jump to content

Basic Graphing Script

From Luter 345 Experiments
Revision as of 15:40, 18 September 2018 by en>Brash (ROOT Script)
(diff) ← Older revision | Latest revision (diff) | Newer revision → (diff)

ROOT Script (must be called graph_basic.C)

//
// Fit to U.S. Standard Atmosphere Data
// From http://www.engineeringtoolbox.com/standard-atmosphere-d_604.html
//

#include <ostream>

Double_t fitfunction(Double_t *x, Double_t *par)
{
        Double_t f = par[0]*exp(par[1]*x[0]+par[2]*x[0]*x[0]);
        return f;
}

void graph_basic(Bool_t ylog = 0, Bool_t xlog = 0, Bool_t ifit = 0) {

   gStyle->SetOptFit(1);
   TCanvas *c1 = new TCanvas("c1","Density of Air vs. Altitude",200,10,1000,750);
   FILE *f = fopen("density.txt","r");

   c1->SetFillColor(42);
   c1->SetGrid();
   if (xlog) {c1->SetLogx();}
   if (ylog) {c1->SetLogy();}

   const Int_t n = 10000;
   Double_t altitude[n], temp[n], gravity[n], pressure[n], density[n], viscosity[n];
   Int_t i=0;

   while (!feof(f)){
          fscanf(f,"%lf %lf %lf %lf %lf %lf\n",&altitude[i],&temp[i],&gravity[i],&pressure[i],&density[i],&viscosity[i]);
          density[i]=density[i]*0.1;
          i++;
   }

   const Int_t np=i;

   TGraph *gr = new TGraph(np,altitude,density);
   TF1 *pfit1 = new TF1("fitfunction",fitfunction,0.0,80000.0,3);
   pfit1->SetParameters(1.2,-0.000000001,.00000000001);
   pfit1->SetParNames("rho_0","linear","quadratic");

   gr->SetMarkerColor(kBlue);
   gr->SetMarkerStyle(21);
   gr->SetTitle("Density");
   gr->GetXaxis()->SetTitle("Altitude (m)");
   gr->GetYaxis()->SetTitle("Density (kg/m^3)");

//   TAxis *axis = gr->GetXaxis();
//   axis->SetLimits(-5.,5.);
//   gr->GetHistogram()->SetMaximum(40000.0);
//   gr->GetHistogram()->SetMinimum(0.0);

   if (ifit) gr->Fit("fitfunction","V");
   gr->Draw("AP");

   // TCanvas::Update() draws the frame, after which one can change it
   c1->Update();
   c1->GetFrame()->SetFillColor(21);
   c1->GetFrame()->SetBorderSize(12);
   c1->Modified();
}

Data File (density.txt)

0       15.00   9.807   10.13   12.25   1.789
1000    8.50    9.804   8.988   11.12   1.758
2000    2.00    9.801   7.950   10.07   1.726
3000    -4.49   9.797   7.012   9.093   1.694
4000    -10.98  9.794   6.166   8.194   1.661
5000    -17.47  9.791   5.405   7.364   1.628
6000    -23.96  9.788   4.722   6.601   1.595
7000    -30.45  9.785   4.111   5.900   1.561
8000    -36.94  9.782   3.565   5.258   1.527
9000    -43.42  9.779   3.080   4.671   1.493
10000   -49.90  9.776   2.650   4.135   1.458
15000   -56.50  9.761   1.211   1.948   1.422
20000   -56.50  9.745   0.5529  0.8891  1.422
25000   -51.60  9.730   0.2549  0.4008  1.448
30000   -46.64  9.715   0.1197  0.1841  1.475
40000   -22.80  9.684   0.0287  0.03996 1.601
50000   -25.00  9.654   0.007978        0.01027 1.704
60000   -26.13  9.624   0.002196        0.003097        1.584
70000   -53.57  9.594   0.00052 0.0008283       1.438
80000   -74.51  9.564   0.00011 0.0001846       1.321