]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add a class which can generate 2d and 3d histograms.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 10 Mar 1999 09:57:11 +0000 (09:57 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 10 Mar 1999 09:57:11 +0000 (09:57 +0000)
git-svn-id: https://svn.dealii.org/trunk@980 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/histogram.h [new file with mode: 0644]
deal.II/deal.II/source/numerics/histogram.cc [new file with mode: 0644]

diff --git a/deal.II/deal.II/include/numerics/histogram.h b/deal.II/deal.II/include/numerics/histogram.h
new file mode 100644 (file)
index 0000000..172976d
--- /dev/null
@@ -0,0 +1,217 @@
+/*----------------------------   histogram.h     ---------------------------*/
+/*      $Id$                 */
+#ifndef __histogram_H
+#define __histogram_H
+/*----------------------------   histogram.h     ---------------------------*/
+
+
+#include <lac/forward-declarations.h>
+#include <vector>
+
+
+
+
+/**
+ * This class provides some facilities to generate 2d and 3d histograms.
+ * It is used by giving it one or several data sets and a rule how to
+ * break the range of values therein into intervals (e.g. linear spacing
+ * or logarithmic spacing of intervals). The values are then sorted into
+ * the different intervals and the number of values in each interval is
+ * stored for output later. In case only one data set was given, the
+ * resulting histogram will be a 2d one, while it will be a 3d one if
+ * more than one data set was given. For more than one data set, the same
+ * intervals are used for each of them anyway, to make comparison easier.
+ *
+ *
+ * \subsection{Ways to generate the intervals}
+ *
+ * At present, the following schemes for interval spacing are implemented:
+ * \begin{itemize}
+ * \item Linear spacing: The intervals are distributed in constant steps
+ *    between the minimum and maximum values of the data.
+ * \item Logaritmic spacing: The intervals are distributed in constant
+ *    steps between the minimum and maximum values of the logs of the values.
+ *    This scheme is only useful if the data has only positive values.
+ *    Negative and zero values are sorted into the leftmost interval.
+ * \end{itemize}
+ *
+ *
+ * \subsection{Output formats}
+ *
+ * At present, only GNUPLOT output is supported.
+ *
+ *
+ * @author Wolfgang Bangerth, 1999
+ */
+class Histogram
+{
+  public:
+                                    /**
+                                     * Definition of several ways to arrange
+                                     * the spacing of intervals.
+                                     */
+    enum IntervalSpacing {
+         linear, logarithmic
+    };
+    
+
+                                    /**
+                                     * Take several lists of values, each on
+                                     * to produce one histogram that will
+                                     * then be arrange one behind each other.
+                                     *
+                                     * Using several data sets at once allows
+                                     * to compare them more easily, since
+                                     * the intervals into which the data is
+                                     * sorted is the same for all data sets.
+                                     *
+                                     * The histograms will be arranged such
+                                     * that the computed intervals of the
+                                     * #values[i][j]# form the x-range,
+                                     * and the number of values in each
+                                     * interval will be the y-range (for
+                                     * 2d plots) or the z-range (for 3d
+                                     * plots). For 3d plots, the #y_values#
+                                     * parameter is used to assign each
+                                     * data set a value in the y direction,
+                                     * which is the depth coordinate in the
+                                     * resulting plot. For 2d plots,
+                                     * the #y_values# are ignored.
+                                     *
+                                     * If you give only one data set, i.e.
+                                     * #values.size()==1#, then the
+                                     * resulting histogram will be a 2d
+                                     * one.
+                                     *
+                                     * #n_intervals# denotes the number of
+                                     * intervals into which the data will be
+                                     * sorted; #interval_spacing# denotes
+                                     * the way the bounds of the intervals
+                                     * are computed. Refer to the general
+                                     * documentation for more information
+                                     * on this.
+                                     */
+    template <typename number>
+    void evaluate (const vector<Vector<number> > &values,
+                  const vector<double>          &y_values,
+                  const unsigned int             n_intervals,
+                  const IntervalSpacing          interval_spacing = linear);
+
+                                    /**
+                                     * This function is only a wrapper
+                                     * to the above one in case you have
+                                     * only one data set.
+                                     */
+    template <typename number>
+    void evaluate (const Vector<number>          &values,
+                  const unsigned int             n_intervals,
+                  const IntervalSpacing          interval_spacing = linear);
+
+                                    /**
+                                     * Write the histogram computed by
+                                     * the #evaluate# function to a
+                                     * stream in a format suitable to
+                                     * the GNUPLOT program. The function
+                                     * generates 2d or 3d histograms.
+                                     */
+    void write_gnuplot (ostream &out) const;
+
+                                    /**
+                                     * Exception.
+                                     */
+    DeclException0 (ExcEmptyData);
+                                    /**
+                                     * Exception.
+                                     */
+    DeclException0 (ExcInvalidIntervals);
+                                    /**
+                                     * Exception.
+                                     */
+    DeclException0 (ExcInvalidData);
+                                    /**
+                                     * Exception.
+                                     */
+    DeclException0 (ExcIO);
+                                    /**
+                                     * Exception.
+                                     */
+    DeclException2 (ExcIncompatibleArraySize,
+                   int, int,
+                   << "The two array sizes " << arg1 << " and " << arg2
+                   << " must match, but don't.");
+    
+  private:
+                                    /**
+                                     * Structure denoting one
+                                     * interval.
+                                     */
+    struct Interval 
+    {
+                                        /**
+                                         * Constructor. Sets the
+                                         * bounds and sets the number of
+                                         * values in this interval to zero.
+                                         */
+       Interval (const double left_point,
+                 const double right_point);
+
+                                        /**
+                                         * Left bound of the interval.
+                                         */
+       double       left_point;
+
+                                        /**
+                                         * Right bound of the interval.
+                                         */
+       double       right_point;
+
+                                        /**
+                                         * Number of values in this
+                                         * interval.
+                                         */
+       unsigned int content;
+    };
+
+                                    /**
+                                     * "Less-than" operation which
+                                     * finds the minimal positive value
+                                     * by sorting zero and negative value
+                                     * to be larger than the largest positive
+                                     * number. Used to find the lower bound
+                                     * of the leftmost interval in the
+                                     * logarithmic case interval spacing
+                                     * scheme.
+                                     *
+                                     * Return #true#, if (#n1<n2#,
+                                     * and (#n1>0# or #n2<0#)), or
+                                     * (n2<n1 and n1>0 and n2<=0). This
+                                     * in effect sorts all negativ
+                                     * numbers to be larger than the
+                                     * largest positive number.
+                                     */
+    template <typename number>
+    static bool logarithmic_less (const number n1,
+                                 const number n2);
+    
+                                    /**
+                                     * Vector holding one set of intervals
+                                     * for each data set given to the
+                                     * #evaluate# function.
+                                     */
+    vector<vector<Interval> > intervals;
+
+                                    /**
+                                     * Values for the depth axis of 3d
+                                     * histograms. Stored in the #evaluate#
+                                     * function.
+                                     */
+    vector<double>            y_values;
+};
+
+
+
+
+/*----------------------------   histogram.h     ---------------------------*/
+/* end of #ifndef __histogram_H */
+#endif
+/*----------------------------   histogram.h     ---------------------------*/
diff --git a/deal.II/deal.II/source/numerics/histogram.cc b/deal.II/deal.II/source/numerics/histogram.cc
new file mode 100644 (file)
index 0000000..030872d
--- /dev/null
@@ -0,0 +1,266 @@
+/* $Id$ */
+
+#include <lac/vector.h>
+#include <basic/histogram.h>
+#include <algorithm>
+#include <cmath>
+
+
+
+
+
+template <typename number>
+inline
+bool Histogram::logarithmic_less (const number n1,
+                                 const number n2)
+{
+  return (((n1<n2) && (n1>0)) ||
+         ((n1<n2) && (n2<=0)) ||
+         ((n2<n1) && (n1>0) && (n2<=0)));
+};
+
+
+
+Histogram::Interval::Interval (const double left_point,
+                              const double right_point) :
+               left_point (left_point),
+               right_point (right_point),
+               content (0)
+{};
+
+
+
+template <typename number>
+void Histogram::evaluate (const vector<Vector<number> > &values,
+                         const vector<double>          &_y_values,
+                         const unsigned int             n_intervals,
+                         const IntervalSpacing          interval_spacing)
+{
+  Assert (values.size() > 0, ExcEmptyData());
+  Assert (n_intervals > 0, ExcInvalidIntervals());
+  for (unsigned int i=0; i<values.size(); ++i)
+    Assert (values[i].size() > 0, ExcEmptyData());
+  Assert (values.size() == _y_values.size(),
+         ExcIncompatibleArraySize(values.size(), _y_values.size()));
+
+                                  // store y_values
+  y_values = _y_values;
+  
+                                  // first find minimum and maximum value
+                                  // in the indicators
+  number min_value, max_value;
+  switch (interval_spacing)
+    {
+      case linear:
+           min_value = *min_element(values[0].begin(),
+                                    values[0].end());
+           max_value = *max_element(values[0].begin(),
+                                    values[0].end());
+
+           for (unsigned int i=1; i<values.size(); ++i)
+             {
+               min_value = min (min_value,
+                                *min_element(values[i].begin(),
+                                             values[i].end()));
+               max_value = max (max_value,
+                                *max_element(values[i].begin(),
+                                             values[i].end()));
+             };
+           
+           break;
+           
+      case logarithmic:
+           min_value = *min_element(values[0].begin(),
+                                    values[0].end(),
+                                    &logarithmic_less<number>);
+           max_value = *max_element(values[0].begin(),
+                                    values[0].end(),
+                                    &logarithmic_less<number>);
+
+           for (unsigned int i=1; i<values.size(); ++i)
+             {
+               min_value = min (min_value,
+                                *min_element(values[i].begin(),
+                                             values[i].end(),
+                                             &logarithmic_less<number>),
+                                &logarithmic_less<number>);
+               max_value = max (max_value,
+                                *max_element(values[i].begin(),
+                                             values[i].end(),
+                                             &logarithmic_less<number>));
+             };
+           
+           break;
+
+      default:
+           Assert (false, ExcInternalError());
+    };
+  
+  Assert (max_value > min_value, ExcInvalidData());
+  
+
+
+                                  // now set up the intervals based on
+                                  // the min and max values
+  intervals.clear ();
+                                  // set up one list of intervals
+                                  // for the first data vector. we will
+                                  // then produce all the other lists
+                                  // for the other data vectors by
+                                  // copying
+  intervals.push_back (vector<Interval>());
+  
+  switch (interval_spacing)
+    {
+      case linear:
+      {
+       const float delta = (max_value-min_value)/n_intervals;
+
+       for (unsigned int n=0; n<n_intervals; ++n)
+         intervals[0].push_back (Interval(min_value+n*delta,
+                                          min_value+(n+1)*delta));
+
+       break;
+      };
+
+      case logarithmic:
+      {
+       const float delta = (log(max_value)-log(min_value))/n_intervals;
+
+       for (unsigned int n=0; n<n_intervals; ++n)
+         intervals[0].push_back (Interval(exp(log(min_value)+n*delta),
+                                          exp(log(min_value)+(n+1)*delta)));
+       
+       break;
+      };
+
+      default:
+           Assert (false, ExcInternalError());
+    };
+
+                                  // fill the other lists of intervals
+  for (unsigned int i=1; i<values.size(); ++i)
+    intervals.push_back (intervals[0]);
+
+
+                                  // finally fill the intervals
+  for (unsigned int i=0; i<values.size(); ++i)
+    for (Vector<number>::const_iterator p=values[i].begin();
+        p < values[i].end(); ++p)
+      {
+                                        // find the right place for *p in
+                                        // intervals[i]. use regular
+                                        // operator< here instead of
+                                        // the logarithmic one to
+                                        // map negative or zero value
+                                        // to the leftmost interval always
+       for (unsigned int n=0; n<n_intervals; ++n)
+         if (*p <= intervals[i][n].right_point)
+           {
+             ++intervals[i][n].content;
+             break;
+           };
+      };
+};
+
+
+
+template <typename number>
+void Histogram::evaluate (const Vector<number>    &values,
+                         const unsigned int       n_intervals,
+                         const IntervalSpacing    interval_spacing) 
+{
+  vector<Vector<number> > values_list (1,
+                                      values);
+  evaluate (values_list, vector<double>(1,0.), n_intervals, interval_spacing);
+};
+
+
+
+void Histogram::write_gnuplot (ostream &out) const
+{
+  AssertThrow (out, ExcIO());
+  Assert (intervals.size() > 0, ExcEmptyData());
+
+                                  // do a simple 2d plot, if only
+                                  // one data set is available
+  if (intervals.size()==1)
+    {
+      for (unsigned int n=0; n<intervals[0].size(); ++n)
+       out << intervals[0][n].left_point
+           << ' '
+           << intervals[0][n].content
+           << endl
+           << intervals[0][n].right_point
+           << ' '
+           << intervals[0][n].content
+           << endl;
+    }
+  else
+                                    // otherwise create a whole 3d plot
+                                    // for the data. use th patch method
+                                    // of gnuplot for this
+    for (unsigned int i=0; i<intervals.size(); ++i)
+      {
+       for (unsigned int n=0; n<intervals[i].size(); ++n)
+         out << intervals[i][n].left_point
+             << ' '
+             << i
+             << ' '
+             << intervals[i][n].content
+             << endl
+             << intervals[i][n].right_point
+             << ' '
+             << i
+             << ' '
+             << intervals[i][n].content
+             << endl;
+       
+       out << endl;
+       
+       for (unsigned int n=0; n<intervals[i].size(); ++n)
+         out << intervals[i][n].left_point
+             << ' '
+             << i+1
+             << ' '
+             << intervals[i][n].content
+             << endl
+             << intervals[i][n].right_point
+             << ' '
+             << i+1
+             << ' '
+             << intervals[i][n].content
+             << endl;
+
+       out << endl;
+      };
+
+  AssertThrow (out, ExcIO());
+};
+
+
+
+
+// explicit instantiations for float
+template
+void Histogram::evaluate (const vector<Vector<float> >  &values,
+                         const vector<double>          &y_values, 
+                         const unsigned int             n_intervals,
+                         const IntervalSpacing          interval_spacing);
+template
+void Histogram::evaluate (const Vector<float>   &values,
+                         const unsigned int     n_intervals,
+                         const IntervalSpacing  interval_spacing);
+
+
+// explicit instantiations for double
+template
+void Histogram::evaluate (const vector<Vector<double> >  &values,
+                         const vector<double>           &y_values, 
+                         const unsigned int              n_intervals,
+                         const IntervalSpacing           interval_spacing);
+template
+void Histogram::evaluate (const Vector<double>   &values,
+                         const unsigned int      n_intervals,
+                         const IntervalSpacing   interval_spacing);
+

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.