graph_histograms.hh 5.98 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
// Copyright (C) 2008  Tiago de Paula Peixoto <tiago@forked.de>
//
// This program is free software; you can redistribute it and/or
// modify it under the terms of the GNU General Public License
// as published by the Free Software Foundation; either version 3
// of the License, or (at your option) any later version.
//
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
// GNU General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with this program. If not, see <http://www.gnu.org/licenses/>.

16
17
18
19
20
21
22
23
24
25
26
27
28
29
#ifndef GRAPH_HISTOGRAMS_HH
#define GRAPH_HISTOGRAMS_HH

#include <algorithm>
#include <boost/numeric/conversion/bounds.hpp>
#include <boost/numeric/conversion/cast.hpp>
#include <boost/python/object.hpp>
#include <boost/python/list.hpp>
#include <boost/python/extract.hpp>
#include "numpy_bind.hh"
#include "histogram.hh"
#include "shared_map.hh"

namespace graph_tool
30
{
31
32
using namespace std;
using namespace boost;
33

34
class VertexHistogramFiller
35
{
36
37
38
39
40
public:
    template <class Graph, class DegreeSelector, class ValueType>
    void UpdateRange(Graph& g,
                     typename graph_traits<Graph>::vertex_descriptor v,
                     DegreeSelector& deg, pair<ValueType,ValueType>& range)
41
    {
42
43
44
        ValueType val = deg(v, g);
        range.first = min(range.first, val);
        range.second = max(range.second, val);
45
    }
46
47
48
49

    template <class Graph, class DegreeSelector, class Hist>
    void operator()(Graph& g, typename graph_traits<Graph>::vertex_descriptor v,
                    DegreeSelector& deg, Hist& hist)
50
    {
51
52
53
        typename Hist::point_t p;
        p[0] = deg(v, g);
        hist.PutValue(p);
54
    }
55
};
56

57
class EdgeHistogramFiller
58
{
59
60
61
62
63
public:
    template <class Graph, class EdgeProperty, class ValueType>
    void UpdateRange(Graph& g,
                     typename graph_traits<Graph>::vertex_descriptor v,
                     EdgeProperty& eprop, pair<ValueType,ValueType>& range)
64
    {
65
66
67
        typename graph_traits<Graph>::out_edge_iterator e, e_begin, e_end;
        tie(e_begin,e_end) = out_edges(v, g);
        for(e = e_begin; e != e_end; ++e)
68
        {
69
70
71
            ValueType val = eprop[*e];
            range.first = min(range.first, val);
            range.second = max(range.second, val);
72
        }
73
    }
74

75
76
77
78
    template <class Graph, class EdgeProperty, class Hist>
    void operator()(Graph& g, typename graph_traits<Graph>::vertex_descriptor v,
                    EdgeProperty& eprop, Hist& hist)
    {
79
        typename Hist::point_t p;
80
81
82
        typename graph_traits<Graph>::out_edge_iterator e, e_begin, e_end;
        tie(e_begin,e_end) = out_edges(v,g);
        for(e = e_begin; e != e_end; ++e)
83
        {
84
85
            p[0] = eprop[*e];
            hist.PutValue(p);
86
87
88
89
        }
    }
};

90
91
92
// generalized functor to obtain histogram of different types of "degrees"
template <class HistogramFiller>
struct get_histogram
93
{
94
95
96
97
98
    get_histogram(python::object& hist, const vector<long double>& bins,
                  python::object& ret_bins)
        : _hist(hist), _bins(bins), _ret_bins(ret_bins) {}

    template <class Graph, class DegreeSelector>
99
    void operator()(Graph& g, DegreeSelector deg) const
100
    {
101
102
        typedef typename DegreeSelector::value_type value_type;
        typedef Histogram<value_type, size_t, 1> hist_t;
103

104
        HistogramFiller filler;
105

Tiago Peixoto's avatar
Tiago Peixoto committed
106
        vector<value_type> bins(_bins.size());
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
        for (size_t i = 0; i < bins.size(); ++i)
        {
            // we'll attempt to recover from out of bounds conditions
            try
            {
                bins[i] = numeric_cast<value_type,long double>(_bins[i]);
            }
            catch (boost::numeric::negative_overflow&)
            {
                bins[i] = boost::numeric::bounds<value_type>::lowest();
            }
            catch (boost::numeric::positive_overflow&)
            {
                bins[i] = boost::numeric::bounds<value_type>::highest();
            }
        }
123

124
125
        // sort the bins
        sort(bins.begin(), bins.end());
126

127
128
129
130
131
132
133
134
135
        // clean bins of zero size
        vector<value_type> temp_bin(1);
        temp_bin[0] = bins[0];
        for (size_t j = 1; j < bins.size(); ++j)
        {
            if (bins[j] > bins[j-1])
                temp_bin.push_back(bins[j]);
        }
        bins = temp_bin;
136

137
138
        // find the data range
        pair<value_type,value_type> range;
139

140
141
142
143
144
145
146
147
148
149
150
151
152
        if (bins.size() == 1)
        {
            typename graph_traits<Graph>::vertex_iterator vi,vi_end;
            range.first = boost::numeric::bounds<value_type>::highest();
            range.second = boost::numeric::bounds<value_type>::lowest();
            for (tie(vi, vi_end) = vertices(g); vi != vi_end; ++vi)
                filler.UpdateRange(g, *vi, deg, range);
        }
        else
        {
            range.first = bins.front();
            range.second = bins.back();
        }
153
154
155
156
157
158
159
160

        boost::array<pair<value_type, value_type>, 1> data_range;
        data_range[0] = range;
        boost::array<vector<value_type>, 1> bin_list;
        bin_list[0] = bins;

        hist_t hist(bin_list, data_range);
        SharedHistogram<hist_t> s_hist(hist);
161
162

        int i, N = num_vertices(g);
163
164
        #pragma omp parallel for default(shared) private(i) \
            firstprivate(s_hist) schedule(dynamic)
165
166
167
168
169
        for (i = 0; i < N; ++i)
        {
            typename graph_traits<Graph>::vertex_descriptor v = vertex(i, g);
            if (v == graph_traits<Graph>::null_vertex())
                continue;
Tiago Peixoto's avatar
Tiago Peixoto committed
170
            filler(g, v, deg, s_hist);
171
        }
172
        s_hist.Gather();
173

174
175
176
177
        bin_list = hist.GetBins();
        python::object ret_bins = wrap_vector_owned(bin_list[0]);
        _ret_bins = ret_bins;
        _hist = wrap_multi_array_owned<size_t,1>(hist.GetArray());
178
    }
179
180
181
182
    python::object& _hist;
    const vector<long double>& _bins;
    python::object& _ret_bins;
};
183

184
} // graph_tool namespace
185

186
#endif // GRAPH_HISTOGRAMS_HH