Jak na grafy v Gnuplotu

?0
Přichází další ze série "quickguide" nebo "tahák", tentokráte pro Gnuplot, protože jsem nedávno potřeboval nějak jednoduše z programu v C++ vykreslit graf a nenašel jsem jednoduchou knihovnu, která by to uměla. Proto jsem spíš než po velké komplexní knihovně sáhl po Gnuplotu, o kterém vím již dlouho, jen jsem doposud neměl touhu se mu věnovat, protože to není zrovna jako klikání grafu v Excelu . Je na něm vidět, že je to nástroj "linuxářů". Nicméně se jedná o multiplatformní aplikaci.
Tento příspěvek má 4 krátké podkapitoly:

1. Základy

Domovská stránka projektu Gnuplot: http://gnuplot.sourceforge.net/.

1.1. Nápověda

Gnuplot obsahuje vestavěnou nápovědu (na Windows je ve formátu komprimované HTML nápovědy - CHM). Nápovědu spustíme příkazem help. Pokud nás zajímá nějaký konkrétní příkaz, lze ho předat jako argument, např. help set nebo dokonce help set xrange.

1.2. Výstupní prostředí / Terminály

Ač to nebylo úmyslem, tak článek hezky navazuje na předchozí, o matematické sazbě. Gnuplot totiž umožňuje výstup do velkého množství formátů. Kromě výstupu do konzole a do okna v grafickém prostředí (X11,Windows,...) je možné ukládat výsledné grafy do souborů ve formátech PNG, SVG, HTML5 Canvas, PostScipt, TeX a spoustu dalších. Např. zmíněný výstup do TeXu je hodně užitečný při sázení dokumentů. Výsledný graf (output.tex) lze velice snadno vložit pomocí TeXového příkazu \input:

\documentclass{article}

\begin{document}

Graf z GNU Plotu (output = latex):

\input{output.tex}

\end{document}
Stejně tak SVG a HTML5 Canvas jsou užitečné pro zpracování na webu, protože narozdíl od obyčejného obrázku jsou tyto dva formáty interaktivní (v grafu se dá zoomovat, zapínat/vypínat grid, jednoduše odečítat body z grafu, atd.). Způsob výstupu se mění pomocí příkazu set terminal, případně ve zkrácené verzi set term. Např.: set term png, set term tex, set term svg, set term wxt, atd.

1.3. Kreslení grafů funkcí

Některé důležité příkazy:

  • set - slouží k nastavení parametrů grafu, k přepínání módů a typu výstupů - co všechno lze ovlivnit, je uvedeno v nápovědě
  • reset - zruší veškerá nastavení vzniklá použitím příkazu set
  • plot - kreslí graf - zpravidla následuje jméno funkce a způsob vykreslení, viz následující ukázky
  • replot - překreslí křivku do předem připraveného prostředí vzniklého použitím příkazu plot; argumentem je zpravidla funkce, kterou chceme vykreslit

Graf funkce sinus:

  1. nastav rozsah osy X (pi je předdefinovaná konstanta Gnuplotu)
  2. zobraz grid a nastav jeho hustotu (dvojice "název bodu" hodnota)
  3. popiš osy X a Y
  4. vykresli funkci sinus

set xrange [-pi:pi]            # we want only one cycle
set xtics ("0" 0, "90" pi/2, "-90" -pi/2, \
           "" pi/4 1, "" -pi/4 1, "" 3*pi/4 1, \
           "" -3*pi/4 1)
set grid
set xlabel "Angle,\n in degrees"
set ylabel "sin(angle)"
plot sin(x)

Více křivek v jednom grafu (sinus, cosinus):

  1. umísti legendu do levého horního rohu a udělej kolem ní rámeček
  2. vykresli sinus - v každém bodu průběhu vykresli značku (v tomto případě identifikovanou číslem 5 = trojúhelník)
  3. vykresli cosinus a vyplň plochu mezi křivkou a osou X (boxes)
Pozn.: při plotu sin(x) jsou použity plné názvy příkazů, zatímco ve vykreslení cos(x) jsou použity jejich zkrácené ekvivalenty (t = title, w = with, lt = linetype).

set key top left
set key box
plot [-pi:pi] sin(x) title "sinusoid" \
              with linespoints pointtype 5, \
              cos(x) t 'cosine' w boxes lt 4

Multiplot (více nezávislých grafů v jednom výstupu):

  1. Nastav rozsah osy X a hranice grafu
  2. Přepni do módu multiplot
  3. Pro každý graf
    1. Nastav velikost dílčího grafu
    2. Nastav pozici dílčího grafu
    3. Vykresli křivku do dílčího grafu
  4. Přepni zpět do módu singleplot
  5. Vyresetuj nastavení do výchozího stavu

set xrange [-pi:pi]

# Uncomment the following to line up the axes
# set lmargin 6

# Gnuplot recommends setting the size and origin before
# going to multiplot mode
# This sets up bounding boxes and may be required on
# some terminals
set size 1,1
set origin 0,0

# Done interactively, this takes gnuplot into multiplot
# mode and brings up a new prompt ("multiplot >"
# instead of "gnuplot >")
set multiplot

# plot the first graph so that it takes a quarter of
# the screen
set size 0.5,0.5
set origin 0,0.5
plot sin(x)

# plot the second graph so that it takes a quarter of
# the screen
set size 0.5,0.5
set origin 0,0
plot 1/sin(x)

# plot the third graph so that it takes a quarter of
# the screen
set size 0.5,0.5
set origin 0.5,0.5
plot cos(x)

# plot the fourth graph so that it takes a quarter of
# the screen
set size 0.5,0.5
set origin 0.5,0
plot 1/cos(x)

# On some terminals, nothing gets plotted until
# this command is issued
unset multiplot

# remove all customization
reset

Jeden graf a dvěma osy s různými měřítky:

  1. Nastav rozsah osy X1
  2. Nastav rozsah osy Y1
  3. Nastav rozsah osy Y2
  4. Nastav rozestupy mřížky na osách X1, Y1 a Y2
  5. Vykresli sin(x) do systému X1Y1 a druhý sin(x) do systému X1Y2
  6. Čekej neomezenou dobu na stisk Enteru (v grafickém rozhraní vyskočí dialogové okno)

set xrange [0:31.4159] #rozsah osy x
set yrange [-1:1]
set y2range [-2:2] #změněný rozsah pro osu y2
set xtics 3.14159*2#chceme mít značku pro každou celou vlnu
set ytics nomirror #vypneme rušivé značky z osy y na ose y2
set y2tics -2,0.5,2#značky na ose y2 s krokem 0.5

plot sin(x), sin(x) axes x1y2 #první sinus pro x1y1

pause -1 "Stiskni enter..." #čeká na enter, vhodné v xsech


1.4. Výpočet parametrů funkcí / "Fitování"

Gnuplot obsahuje i poměrně pokročilý matematický aparát. Velmi užitečná je funkce fit. Její použití je nejlépe vidět na jednoduchém příkladu:

f(x) = a * sin(b*x) + c           # funkce, kterou fituju
fit f(x) "mereni.dat" via a, b, c # pomoci parametru a,b,c
Výstup vypadá takto:
  After 14 iterations the fit converged.
  final sum of squares of residuals : 8.24812e-31
  abs. change during last iteration : -1.08701e-28


  Hmmmm.... Sum of squared residuals is zero. Can't compute errors.

  Final set of parameters
  =======================

  a               = 1
  b               = 1.5708
  c               = -2.04925e-21
Hodnoty parametrů už známe, takže můžeme zobrazit jak průběh funkce souhlasí, či nesouhlasí s naměřenými hodnotami:
plot "mereni.dat", f(x)

K odhadu parametrů se používá kritérium nejmenších čtverců (Least Squares Error), nad kterým operuje známý algoritmus Levenberg-Marquardt.

2. Příklady grafů


Následuje seznam pokročilých a zároveň pěkných grafů, které budu v brzké budoucnosti potřebovat použít (alespoň si to myslím ). Jedná se o vybrané grafy z ukázek přímo v dokumentaci Gnuplotu (Demo v4.6). Pouze příklad s histogramem je převzat z článku Statistic Analysis and Histogram.
Gnuplot defaultně obsahuje i pár datových souborů, které byly v následujících příkladech použity, konkrétně jde o datasety silver.dat a immigration.dat.

2.1. Heatmap


2.1.1. Flat

#
# Generating a 2D heat map from ascii data
#

set title "Heat Map generated from a file containing Z values only"
unset key
set tic scale 0

# Color runs from white to green
set palette rgbformula -7,2,-7
set cbrange [0:5]
set cblabel "Score"
unset cbtics

set xrange [-0.5:4.5]
set yrange [-0.5:4.5]

set view map
splot '-' matrix with image
5 4 3 1 0
2 2 0 0 1
0 0 0 1 0
0 0 0 2 3
0 1 2 4 3
e
e

2.1.2. 3D

#
# Demonstrate use of 4th data column to color a 3D surface.
# Also demonstrate use of the pseudodata special file '++'.
# This plot is nice for exploring the effect of
# the 'l' and 'L' hotkeys.
#
set view 49, 28, 1, 1.48
set xrange [ 5 : 35 ] noreverse nowriteback
set yrange [ 5 : 35 ] noreverse nowriteback
# set zrange [ 1.0 : 3.0 ] noreverse nowriteback
set ticslevel 0
set format cb "%4.1f"
set colorbox user size .03, .6 noborder
set cbtics scale 0

set samples 25, 25
set isosamples 50, 50

set title "4D data (3D Heat Map)"\
          ."\nIndependent value color-mapped onto 3D surface"\
          offset 0,1
set xlabel "x" offset 3, 0, 0 
set ylabel "y" offset -5, 0, 0
set zlabel "z" offset 2, 0, 0 
set pm3d implicit at s

Z(x,y) = 100. * (sinc(x,y) + 1.5)
sinc(x,y) = sin(sqrt((x-20.)**2+(y-20.)**2))/\
            sqrt((x-20.)**2+(y-20.)**2)
color(x,y) = 10. * (1.1 + sin((x-20.)/5.)*cos((y-20.)/10.))

splot '++' using 1:2:(Z($1,$2)):(color($1,$2))\
           with pm3d title "4 data columns x/y/z/color"

2.2. Scatter


2.2.1. 2D

# set terminal svg size 600,400 dynamic enhanced fname 'arial'
# fsize 10 mousing name "random_1" butt solid 
# set output 'random.1.svg'
set dummy t,y
set format x "%3.2f"
set format y "%3.2f"
set format z "%3.2f"
unset key
set parametric
set samples 1000, 1000
set style function dots
set title "Lattice test for random numbers" 
set xlabel "rand(n) ->" 
set xrange [ 0.00000 : 1.00000 ] noreverse nowriteback
set ylabel "rand(n + 1) ->" 
set yrange [ 0.00000 : 1.00000 ] noreverse nowriteback
set zlabel "rand(n + 2) ->" 
set zrange [ 0.00000 : 1.00000 ] noreverse nowriteback
plot rand(0), rand(0)

2.2.2. 3D

# set terminal svg size 600,400 dynamic enhanced fname 'arial'
# fsize 10 mousing name "random_2" butt solid 
# set output 'random.2.svg'
set dummy u,v
set format x "%3.2f"
set format y "%3.2f"
set format z "%3.2f"
unset key
set parametric
set samples 50, 50
set style function dots
set title "Lattice test for random numbers" 
set xlabel "rand(n) ->" 
set xrange [ 0.00000 : 1.00000 ] noreverse nowriteback
set ylabel "rand(n + 1) ->" 
set yrange [ 0.00000 : 1.00000 ] noreverse nowriteback
set zlabel "rand(n + 2) ->" 
set zrange [ 0.00000 : 1.00000 ] noreverse nowriteback
splot rand(0), rand(0), rand(0)

2.2.3. 2D + 3D

# set terminal svg size 600,400 dynamic enhanced fname 'arial'
# fsize 10 mousing name "random_3" butt solid 
# set output 'random.3.svg'
set dummy u,v
set format x "%3.2f"
set format y "%3.2f"
set format z "%3.2f"
unset key
set parametric
set view 68, 28, 1, 1
set samples 50, 50
set isosamples 30, 30
set contour base
unset clabel
set hidden3d back offset 1 trianglepattern 3 undefined 1\
    altdiagonal bentover
set cntrparam levels discrete 0.1
set style function dots
set ticslevel 0
set ztics 0.00000,0.05 norangelimit
set title "50 random samples from a 2D Gaussian PDF with\n\
    unit variance, zero mean and no dependence" 
set urange [ -3.00000 : 3.00000 ] noreverse nowriteback
set vrange [ -3.00000 : 3.00000 ] noreverse nowriteback
set xrange [ -3.00000 : 3.00000 ] noreverse nowriteback
set yrange [ -3.00000 : 3.00000 ] noreverse nowriteback
set zrange [ -0.200000 : 0.200000 ] noreverse nowriteback
tstring(n) = sprintf("%d random samples from a 2D Gaussian PDF\
             with\nunit variance, zero mean and no dependence", n)
nsamp = 50
i = 50
GPFUN_tstring = "tstring(n) = sprintf(\"%d random samples from a\
                2D Gaussian PDF with\\nunit variance, zero mean\
                and no dependence\", n)"
splot u,v,( 1/(2*pi) * exp(-0.5 * (u**2 + v**2)) )\
      with line lc rgb "black",\
      "random.tmp" using 1:2:(-0.2)\
      with points pointtype 7 lc rgb "black"

2.2.4. 3D + Histogram

print ""
print "Another Monte Carlo simulation"
print ""
print "This is similar to the previous simulation but uses "
print "multivariate zero mean, unit variance normal data by "
print "computing the distance each point is from the origin. "
print "That distribution is known to fit the Maxwell "
print "probability law, as shown."
print ""
reset
nsamp = 3000
# Generate N random data points.
set print "random.tmp"
do for [i=1:nsamp] {
    print sprintf("%8.5g %8.5g %8.5g", invnorm(rand(0)), \
          invnorm(rand(0)), invnorm(rand(0)))
}
unset print
#
set multiplot layout 1,2
#
unset key
rlow = -4.0
rhigh = 4.0
set parametric
set xrange [rlow:rhigh]; set yrange [rlow:rhigh];
set zrange [rlow:rhigh]
set xtics axis nomirror; set ytics axis nomirror;
set ztics axis nomirror;
set border 0
set xyplane at 0
set xzeroaxis lt -1
set yzeroaxis lt -1
set zzeroaxis lt -1
set view 68, 28, 1.4, 0.9
tstring(n) = sprintf("Gaussian 3D cloud of %d random samples\n", n)
set title tstring(nsamp) offset graph 0.15, graph -0.33
splot "random.tmp" every :::::0 with dots

#
fourinvsqrtpi = 2.25675833419103
maxwell(x,a)=a<=0?1/0:x<0?0.0:fourinvsqrtpi*a**3*x*x*exp(-a*a*x*x)
#
unset parametric
unset xzeroaxis; unset yzeroaxis;
set border
set grid
set samples 200
set size 0.47,0.72
set origin 0.44,0.18
tstring(n) = sprintf("Histogram of distance from origin of\n\
             %d multivariate unit variance samples", n)
set title tstring(nsamp) offset graph 0, graph 0.15
set key bmargin right vertical
xlow = 0.0
xhigh = 4.5
binwidth = 20
scale = (binwidth/(xhigh-xlow))
set xrange [0:xhigh]
set yrange [0:0.65]
bin(x) = (1.0/scale)*floor(x*scale)
plot "random.tmp" using (bin(sqrt($1**2+$2**2+$3**2))):\
     (1.0*scale/nsamp) every :::::0 smooth frequency with steps \
           title "scaled bin frequency", \
           maxwell(x, 1/sqrt(2)) with lines title "Maxwell p.d.f."
#

2.3. Boxplot

#
# Boxplot demo
#
reset

print "*** Boxplot demo ***"

set style fill solid 0.25 border -1
set style boxplot outliers pointtype 7
set style data boxplot
set boxwidth  0.5
set pointsize 0.5

unset key
set border 2
set xtics ("A" 1, "B" 2) scale 0.0
set xtics nomirror
set ytics nomirror
set yrange [0:100]

plot 'silver.dat' using (1):2, '' using (2):(5*$3)

2.4. Histogram

n=100	#number of intervals
max=3.	#max value
min=-3.	#min value
width=(max-min)/n	#interval width
#function used to map a value to the intervals
hist(x,width)=width*floor(x/width)+width/2.0
set term png	#output terminal and file
set output "histogram.png"
set xrange [min:max]
set yrange [0:]
#to put an empty boundary around the
#data inside an autoscaled graph.
set offset graph 0.05,0.05,0.05,0.0
set xtics min,(max-min)/5,max
set boxwidth width*0.9
set style fill solid 0.5	#fillstyle
set tics out nomirror
set xlabel "x"
set ylabel "Frequency"
#count and plot
plot "data.dat" u (hist($1,width)):(1.0)\
     smooth freq w boxes lc rgb"green" notitle
Použitý dataset: data.dat.

2.5. "Error plot"

Jedná se o proložení bodů křivkou a zobrazení rozsahu chyb.

set logscale y
plot "silver.dat" t "rate" w errorb, \
               "" smooth sbezier t "bezier"

2.6. "Klasické" grafy


2.6.1. Lomená čára

#
# Example of using histogram modes
#
reset
set title "US immigration from Europe by decade"
set datafile missing "-"
set xtics nomirror rotate by -45 font ",8"
set key noenhanced
#
# First plot using linespoints
set style data linespoints
plot 'immigration.dat' using 2:xtic(1) title columnheader(2), \
for [i=3:22] '' using i title columnheader(i)
#

2.6.2. Sloupcový graf

#
set title "US immigration from Northern Europe\nPlot selected\
           data columns as histogram of clustered boxes"
set auto x
set yrange [0:300000]
set style data histogram
set style histogram cluster gap 1
set style fill solid border -1
set boxwidth 0.9
set xtic rotate by -45 scale 0 font ",8"
#set bmargin 10 
plot 'immigration.dat' using 6:xtic(1) ti col, '' u 12 ti col,\
     '' u 13 ti col, '' u 14 ti col
#

2.6.3 Sloupcový graf 2 - Stacked Histogram

#
# Stacked histograms
#
set title "US immigration from Europe by decade\n\
           Plot as stacked histogram"
set key invert reverse Left outside
set key autotitle columnheader
set yrange [0:7e6]
set auto x
unset xtics
set xtics nomirror rotate by -45 scale 0 font ",8"
set style data histogram
set style histogram rowstacked
set style fill solid border -1
set boxwidth 0.75
#
plot 'immigration.dat' using 2:xtic(1), for [i=3:22] '' using i
#

3. C++ API pro Gnuplot

Jak už jsem psal v úvodu, hledal jsem něco použitelného pro C/C++. Našel jsem API pro C99 i pro C++03. Pro další práci jsem si vybral právě C++ verzi, protože obě byly velmi podobné, pouze se lišily rozhraním.

Zde je ukázka práce s knihovnou:

#include 
#include 
#include 
#include 
#include "gnuplotpp.h"

using namespace std;
using namespace gnuplotpp;

int main( void )
{
 singleplot g;

 g.set("term png" );
 g.set() << gpset( "output", "\"demo1.png\"" ) << endl;

 g.x.label( "x axis" ).format( "%03.2f" ).tics(0.2).mtics(.05);
 g.y.label( "y axis" );

 g.label( "simgle plot example", coord( "graph 0.5", "graph 0.5" ) );
 g.arrow( coord( "graph 0.5", "graph 0.5" ), coord( "graph 1", "graph 1" ) );
 g.plotfile( "data.dat", 1, 2 )<< " t 'data'" << " pt 1" ;
 g.plotfunc( "x" ) ;

 g.save( "demo1.gp" ).exec();

}
Při bližším prozkoumání jsem zjistil, že se jedná pouze o výstup do std::ostream. Výsledkem je tedy dávka pro Gnuplot, stejně jako kdybych ji psal přímo do souboru. Tato dávka je uložena voláním metody gnuplotpp::plot::save. Metoda gnuplotpp::plot::exec pak nedělá nic jiného než, že spustí Gnuplot s parametrem, který obsahuje cestu k vytvořenému dávkovému souboru - je tedy nutné mít binárku Gnuplotu v systémových cestách.

Líbil se Vám článek?

Share |

comments powered by Disqus