Logo ROOT  
Reference Guide
 
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
Loading...
Searching...
No Matches
gr015_smooth.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_graphs
3/// \notebook
4/// \preview Show scatter plot smoothers: ksmooth, lowess, supsmu,
5/// as described in:
6/// Modern Applied Statistics with S-Plus, 3rd Edition
7/// W.N. Venables and B.D. Ripley
8/// Chapter 9: Smooth Regression, Figure 9.1
9///
10/// Example is a set of data on 133 observations of acceleration against time
11/// for a simulated motorcycle accident, taken from Silverman (1985). The data
12/// are read from motorcycle.dat (2 columns of floats)
13///
14/// \macro_image
15/// \macro_code
16/// \author Christian Stratowa, Vienna, Austria
17
18#include "TString.h"
19#include "TInterpreter.h"
20#include <fstream>
21#include "TH1.h"
22#include "TGraphSmooth.h"
23#include "TCanvas.h"
24#include "TSystem.h"
25
26
29
30void DrawSmooth(Int_t pad, const char *title, const char *xt, const char *yt)
31{
32 vC1->cd(pad);
33 TH1F *vFrame = gPad->DrawFrame(0,-130,60,70);
34 vFrame->SetTitle(title);
35 vFrame->SetTitleSize(0.2);
36 vFrame->SetXTitle(xt);
37 vFrame->SetYTitle(yt);
38 grin->Draw("P");
39 grout->DrawClone("LPX");
40}
41
42void gr015_smooth()
43{
44// data taken from R library MASS: mcycle.txt
45 TString dir = gROOT->GetTutorialDir();
46 dir.Append("/visualisation/graphs/");
47 dir.ReplaceAll("/./","/");
48
49// read file and add to fit object
50 Double_t *x = new Double_t[133];
51 Double_t *y = new Double_t[133];
52 Double_t vX, vY;
53 Int_t vNData = 0;
54 ifstream vInput;
55 vInput.open(Form("%smotorcycle.dat",dir.Data()));
56 while (1) {
57 vInput >> vX >> vY;
58 if (!vInput.good()) break;
59 x[vNData] = vX;
60 y[vNData] = vY;
61 vNData++;
62 }//while
63 vInput.close();
64 grin = new TGraph(vNData,x,y);
65
66// draw graph
67 vC1 = new TCanvas("vC1","Smooth Regression",200,10,900,700);
68 vC1->Divide(2,3);
69
70// Kernel Smoother
71// create new kernel smoother and smooth data with bandwidth = 2.0
72 TGraphSmooth *gs = new TGraphSmooth("normal");
73 grout = gs->SmoothKern(grin,"normal",2.0);
74 DrawSmooth(1,"Kernel Smoother: bandwidth = 2.0","times","accel");
75
76// redraw ksmooth with bandwidth = 5.0
77 grout = gs->SmoothKern(grin,"normal",5.0);
78 DrawSmooth(2,"Kernel Smoother: bandwidth = 5.0","","");
79
80// Lowess Smoother
81// create new lowess smoother and smooth data with fraction f = 2/3
82 grout = gs->SmoothLowess(grin,"",0.67);
83 DrawSmooth(3,"Lowess: f = 2/3","","");
84
85// redraw lowess with fraction f = 0.2
86 grout = gs->SmoothLowess(grin,"",0.2);
87 DrawSmooth(4,"Lowess: f = 0.2","","");
88
89// Super Smoother
90// create new super smoother and smooth data with default bass = 0 and span = 0
91 grout = gs->SmoothSuper(grin,"",0,0);
92 DrawSmooth(5,"Super Smoother: bass = 0","","");
93
94// redraw supsmu with bass = 3 (smoother curve)
95 grout = gs->SmoothSuper(grin,"",3);
96 DrawSmooth(6,"Super Smoother: bass = 3","","");
97
98// cleanup
99 delete [] x;
100 delete [] y;
101 delete gs;
102}
int Int_t
Definition RtypesCore.h:45
double Double_t
Definition RtypesCore.h:59
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
#define gROOT
Definition TROOT.h:406
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition TString.cxx:2489
#define gPad
The Canvas class.
Definition TCanvas.h:23
A helper class to smooth TGraph.
A TGraph is an object made of two arrays X and Y with npoints each.
Definition TGraph.h:41
1-D histogram with a float per channel (see TH1 documentation)
Definition TH1.h:645
Basic string class.
Definition TString.h:139
const char * Data() const
Definition TString.h:376
TString & ReplaceAll(const TString &s1, const TString &s2)
Definition TString.h:704
TString & Append(const char *cs)
Definition TString.h:572
Double_t y[n]
Definition legend1.C:17
Double_t x[n]
Definition legend1.C:17