1 /* PSPP - a program for statistical analysis.
2 Copyright (C) 2004, 2008, 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/>. */
19 #include <output/charts/np-plot.h>
21 #include <gsl/gsl_cdf.h>
23 #include <data/casereader.h>
24 #include <data/casewriter.h>
25 #include <libpspp/message.h>
27 #include <output/chart-provider.h>
28 #include <output/charts/cartesian.h>
29 #include <output/charts/plot-chart.h>
31 #include "gl/minmax.h"
34 #define _(msgid) gettext (msgid)
36 /* An NP or DNP plot. */
41 struct casereader *data;
43 /* Copied directly from struct np. */
45 double dns_min, dns_max;
48 double slope, intercept;
49 double y_first, y_last;
50 double x_lower, x_upper;
54 static const struct chart_class np_plot_chart_class;
55 static const struct chart_class dnp_plot_chart_class;
58 make_np_plot (const struct chart_class *class,
59 const struct np *np, const struct casereader *reader,
62 struct np_plot_chart *npp;
67 npp = xmalloc (sizeof *npp);
68 chart_init (&npp->chart, class);
69 npp->label = xstrdup (label);
70 npp->data = casereader_clone (reader);
71 npp->y_min = np->y_min;
72 npp->y_max = np->y_max;
73 npp->dns_min = np->dns_min;
74 npp->dns_max = np->dns_max;
76 /* Slope and intercept of the ideal normal probability line. */
77 npp->slope = 1.0 / np->stddev;
78 npp->intercept = -np->mean / np->stddev;
80 npp->y_first = gsl_cdf_ugaussian_Pinv (1 / (np->n + 1));
81 npp->y_last = gsl_cdf_ugaussian_Pinv (np->n / (np->n + 1));
83 /* Need to make sure that both the scatter plot and the ideal fit into the
85 npp->x_lower = MIN (np->y_min, (npp->y_first - npp->intercept) / npp->slope);
86 npp->x_upper = MAX (np->y_max, (npp->y_last - npp->intercept) / npp->slope);
87 npp->slack = (npp->x_upper - npp->x_lower) * 0.05;
92 /* Creates and returns a normal probability plot corresponding to
93 the calculations in NP and the data in READER, and label the
94 plot with LABEL. The data in READER must have Y-values in
95 value index NP_IDX_Y and NS-values in value index NP_IDX_NS.
97 Returns a null pointer if the data set is empty.
99 The caller retains ownership of NP and READER. */
101 np_plot_create (const struct np *np, const struct casereader *reader,
104 return make_np_plot (&np_plot_chart_class, np, reader, label);
107 /* Creates and returns a detrended normal probability plot
108 corresponding to the calculations in NP and the data in
109 READER, and label the plot with LABEL. The data in READER
110 must have Y-values in value index NP_IDX_Y and DNS-values in
111 value index NP_IDX_DNS.
113 Returns a null pointer if the data set is empty.
115 The caller retains ownership of NP and READER. */
117 dnp_plot_create (const struct np *np, const struct casereader *reader,
120 return make_np_plot (&dnp_plot_chart_class, np, reader, label);
124 np_plot_chart_draw (const struct chart *chart, plPlotter *lp,
125 struct chart_geometry *geom)
127 const struct np_plot_chart *npp = (struct np_plot_chart *) chart;
128 struct casereader *data;
131 chart_write_title (lp, geom, _("Normal Q-Q Plot of %s"), npp->label);
132 chart_write_xlabel (lp, geom, _("Observed Value"));
133 chart_write_ylabel (lp, geom, _("Expected Normal"));
134 chart_write_xscale (lp, geom,
135 npp->x_lower - npp->slack,
136 npp->x_upper + npp->slack, 5);
137 chart_write_yscale (lp, geom, npp->y_first, npp->y_last, 5);
139 data = casereader_clone (npp->data);
140 for (; (c = casereader_read (data)) != NULL; case_unref (c))
141 chart_datum (lp, geom, 0,
142 case_data_idx (c, NP_IDX_Y)->f,
143 case_data_idx (c, NP_IDX_NS)->f);
144 casereader_destroy (data);
146 chart_line (lp, geom, npp->slope, npp->intercept,
147 npp->y_first, npp->y_last, CHART_DIM_Y);
151 dnp_plot_chart_draw (const struct chart *chart, plPlotter *lp,
152 struct chart_geometry *geom)
154 const struct np_plot_chart *dnpp = (struct np_plot_chart *) chart;
155 struct casereader *data;
158 chart_write_title (lp, geom, _("Detrended Normal Q-Q Plot of %s"),
160 chart_write_xlabel (lp, geom, _("Observed Value"));
161 chart_write_ylabel (lp, geom, _("Dev from Normal"));
162 chart_write_xscale (lp, geom, dnpp->y_min, dnpp->y_max, 5);
163 chart_write_yscale (lp, geom, dnpp->dns_min, dnpp->dns_max, 5);
165 data = casereader_clone (dnpp->data);
166 for (; (c = casereader_read (data)) != NULL; case_unref (c))
167 chart_datum (lp, geom, 0, case_data_idx (c, NP_IDX_Y)->f,
168 case_data_idx (c, NP_IDX_DNS)->f);
169 casereader_destroy (data);
171 chart_line (lp, geom, 0, 0, dnpp->y_min, dnpp->y_max, CHART_DIM_X);
175 np_plot_chart_destroy (struct chart *chart)
177 struct np_plot_chart *npp = (struct np_plot_chart *) chart;
179 casereader_destroy (npp->data);
184 static const struct chart_class np_plot_chart_class =
187 np_plot_chart_destroy
190 static const struct chart_class dnp_plot_chart_class =
193 np_plot_chart_destroy