1 /* PSPP - a program for statistical analysis.
2 Copyright (C) 2004, 2009 Free Software Foundation, Inc.
4 This program is free software: you can redistribute it and/or modify
5 it under the terms of the GNU General Public License as published by
6 the Free Software Foundation, either version 3 of the License, or
7 (at your option) any later version.
9 This program is distributed in the hope that it will be useful,
10 but WITHOUT ANY WARRANTY; without even the implied warranty of
11 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 GNU General Public License for more details.
14 You should have received a copy of the GNU General Public License
15 along with this program. If not, see <http://www.gnu.org/licenses/>. */
22 #include <gsl/gsl_histogram.h>
23 #include <gsl/gsl_randist.h>
26 #include <output/charts/plot-hist.h>
27 #include <output/charts/plot-chart.h>
28 #include <output/chart-provider.h>
30 #include <data/variable.h>
31 #include <libpspp/hash.h>
32 #include <output/chart.h>
33 #include <math/histogram.h>
34 #include <math/moments.h>
37 #define _(msgid) gettext (msgid)
39 static const struct chart_class histogram_chart_class;
41 /* Write the legend of the chart */
43 histogram_write_legend (cairo_t *cr, const struct chart_geometry *geom,
44 double n, double mean, double stddev)
46 double y = geom->data_bottom;
51 char *buf = xasprintf ("N = %.2f", n);
52 cairo_move_to (cr, geom->legend_left, y);
53 chart_label (cr, 'l', 'b', geom->font_size, buf);
54 y += geom->font_size * 1.5;
60 char *buf = xasprintf ("Mean = %.1f", mean);
61 cairo_move_to (cr,geom->legend_left, y);
62 chart_label (cr, 'l', 'b', geom->font_size, buf);
63 y += geom->font_size * 1.5;
69 char *buf = xasprintf ("Std. Dev = %.2f", stddev);
70 cairo_move_to (cr, geom->legend_left, y);
71 chart_label (cr, 'l', 'b', geom->font_size, buf);
79 hist_draw_bar (cairo_t *cr, const struct chart_geometry *geom,
80 const gsl_histogram *h, int bar)
86 const size_t bins = gsl_histogram_bins (h);
87 const double x_pos = (geom->data_right - geom->data_left) * bar / (double) bins ;
88 const double width = (geom->data_right - geom->data_left) / (double) bins ;
90 assert ( 0 == gsl_histogram_get_range (h, bar, &lower, &upper));
92 assert ( upper >= lower);
94 height = gsl_histogram_get (h, bar) *
95 (geom->data_top - geom->data_bottom) / gsl_histogram_max_val (h);
97 cairo_rectangle (cr, geom->data_left + x_pos, geom->data_bottom,
100 cairo_set_source_rgb (cr,
101 geom->fill_colour.red / 255.0,
102 geom->fill_colour.green / 255.0,
103 geom->fill_colour.blue / 255.0);
104 cairo_fill_preserve (cr);
108 draw_tick (cr, geom, TICK_ABSCISSA,
109 x_pos + width / 2.0, "%g", (upper + lower) / 2.0);
112 struct histogram_chart
115 gsl_histogram *gsl_hist;
123 /* Plots a histogram of the data in HIST with the given LABEL.
124 Labels the histogram with each of N, MEAN, and STDDEV that is
125 not SYSMIS. If all three are not SYSMIS and SHOW_NORMAL is
126 true, also draws a normal curve on the histogram. */
128 histogram_chart_create (const struct histogram *hist, const char *label,
129 double n, double mean, double stddev,
132 struct histogram_chart *h;
134 h = xmalloc (sizeof *h);
135 chart_init (&h->chart, &histogram_chart_class);
136 h->gsl_hist = hist->gsl_hist ? gsl_histogram_clone (hist->gsl_hist) : NULL;
137 h->label = xstrdup (label);
141 h->show_normal = show_normal;
146 histogram_chart_draw (const struct chart *chart, cairo_t *cr,
147 struct chart_geometry *geom)
149 struct histogram_chart *h = (struct histogram_chart *) chart;
153 chart_write_title (cr, geom, _("HISTOGRAM"));
155 chart_write_ylabel (cr, geom, _("Frequency"));
156 chart_write_xlabel (cr, geom, h->label);
158 if (h->gsl_hist == NULL)
160 /* Probably all values are SYSMIS. */
164 bins = gsl_histogram_bins (h->gsl_hist);
166 chart_write_yscale (cr, geom, 0, gsl_histogram_max_val (h->gsl_hist), 5);
168 for (i = 0; i < bins; i++)
169 hist_draw_bar (cr, geom, h->gsl_hist, i);
171 histogram_write_legend (cr, geom, h->n, h->mean, h->stddev);
174 && h->n != SYSMIS && h->mean != SYSMIS && h->stddev != SYSMIS)
176 /* Draw the normal curve */
178 double x_min, x_max, not_used;
179 double abscissa_scale;
180 double ordinate_scale;
183 gsl_histogram_get_range (h->gsl_hist, 0, &x_min, ¬_used);
184 range = not_used - x_min;
185 gsl_histogram_get_range (h->gsl_hist, bins - 1, ¬_used, &x_max);
187 abscissa_scale = (geom->data_right - geom->data_left) / (x_max - x_min);
188 ordinate_scale = (geom->data_top - geom->data_bottom) /
189 gsl_histogram_max_val (h->gsl_hist);
191 cairo_move_to (cr, geom->data_left, geom->data_bottom);
192 for (d = geom->data_left;
193 d <= geom->data_right;
194 d += (geom->data_right - geom->data_left) / 100.0)
196 const double x = (d - geom->data_left) / abscissa_scale + x_min;
197 const double y = h->n * range *
198 gsl_ran_gaussian_pdf (x - h->mean, h->stddev);
200 cairo_line_to (cr, d, geom->data_bottom + y * ordinate_scale);
209 histogram_chart_destroy (struct chart *chart)
211 struct histogram_chart *h = (struct histogram_chart *) chart;
212 if (h->gsl_hist != NULL)
213 gsl_histogram_free (h->gsl_hist);
218 static const struct chart_class histogram_chart_class =
220 histogram_chart_draw,
221 histogram_chart_destroy