graph_tree_cts.cc 7.67 KB
Newer Older
1
2
// graph-tool -- a general graph modification and manipulation thingy
//
Tiago Peixoto's avatar
Tiago Peixoto committed
3
// Copyright (C) 2006-2015 Tiago de Paula Peixoto <tiago@skewed.de>
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
//
// 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/>.

#include "graph.hh"

#include "graph_filtering.hh"

#include "graph_selectors.hh"
#include "graph_properties.hh"

Tiago Peixoto's avatar
Tiago Peixoto committed
25
#include <boost/mpl/quote.hpp>
26
#include <boost/graph/breadth_first_search.hpp>
Tiago Peixoto's avatar
Tiago Peixoto committed
27

28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
#include <cmath>

using namespace std;
using namespace boost;
using namespace graph_tool;

typedef pair<double, double> point_t;

point_t interpolate(const point_t& p1, const point_t& p2, double r = 0.5)
{
    point_t ret;
    ret.first = (1 - r) * p1.first + p2.first * r;
    ret.second = (1 - r) * p1.second + p2.second * r;
    return ret;
}

44
45
46
47
48
double dist(point_t& p1, point_t& p2)
{
    return sqrt(pow(p1.first - p2.first, 2) + pow(p1.second - p2.second, 2));
}

49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
void to_bezier(const vector<point_t> &x, vector<point_t>& ncp)
{
    vector<point_t> cp(x.size() + 6);
    for (size_t i = 0; i < 3; ++i)
        cp[i] = x[0];
    for (size_t i = 0; i < x.size(); ++i)
        cp[i + 3] = x[i];
    for (size_t i = cp.size() - 3; i < cp.size(); ++i)
        cp[i] = x.back();

    vector<point_t> one_thirds(cp.size() - 1);
    vector<point_t> two_thirds(cp.size() - 1);

    for (size_t i = 0; i < cp.size() - 1; ++i)
    {
        const point_t& p1 = cp[i];
        const point_t& p2 = cp[i + 1];
        one_thirds[i] = interpolate(p1, p2, 1./3);
        two_thirds[i] = interpolate(p2, p1, 1./3);
    }

    ncp.resize((cp.size() - 3) * 3);
    for (size_t i = 0; i < cp.size() - 3; ++i)
    {
        size_t pos = i * 3;
        ncp[pos] = one_thirds[i + 1];
        ncp[pos + 1] = two_thirds[i + 1];
        ncp[pos + 2] = interpolate(two_thirds[i + 1], one_thirds[i + 2]);
    }
}

void transform(vector<point_t>& cp)
{
    point_t origin = cp[0];
83
    for (auto& xy : cp)
84
    {
85
86
        xy.first -= origin.first;
        xy.second -= origin.second;
87
88
89
90
91
    }

    double t = atan2(cp.back().second - cp.front().second,
                     cp.back().first - cp.front().first);

92
    for (auto& xy : cp)
93
    {
94
95
        double x = xy.first;
        double y = xy.second;
96

97
98
        xy.first = cos(t) * x + sin(t) * y;
        xy.second = -sin(t) * x + cos(t) * y;
99
100
101
102
103
104
105
    }

    point_t d;
    d.first = cp.back().first - cp.front().first;
    d.second = cp.back().second - cp.front().second;
    double r = sqrt(d.first * d.first + d.second * d.second);

106
107
    for (auto& xy : cp)
        xy.first /= r;
108
109
110

    d.first = d.second = 0;
    cp.insert(cp.begin(), d);
111
112
113
}

template <class PosProp>
114
void get_control_points(vector<size_t>& path, PosProp& pos, double beta,
115
116
117
118
119
                        vector<point_t>& ncp)
{
    size_t L = path.size();
    vector<point_t> cp(L);
    for (size_t i = 0; i < L; ++i)
120
121
122
123
124
125
126
127
    {
        auto& p = pos[path[i]];
        if (p.size() < 2)
            p.resize(2);
        cp[i].first = p[0];
        cp[i].second = p[1];
    }

128
129
130
131
132
133
134
135
136
    ncp.resize(L);
    for (size_t i = 0; i < L; ++i)
    {
        ncp[i].first = beta * cp[i].first + (1 - beta) * (cp[0].first + (cp.back().first - cp[0].first) * i / (L - 1.));
        ncp[i].second = beta * cp[i].second + (1 - beta) * (cp[0].second + (cp.back().second - cp[0].second) * i / (L - 1.));
    }
}

template <class Graph>
137
138
void tree_path(Graph& g, size_t s, size_t t, vector<size_t>& path,
               size_t max_depth)
139
140
141
142
143
144
145
146
147
{
    vector<size_t> s_root;
    vector<size_t> t_root;
    s_root.push_back(s);
    t_root.push_back(t);

    size_t v = s;
    size_t u = t;

148
    while (v != u && s_root.size() < max_depth)
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
    {
        typename graph_traits<Graph>::in_edge_iterator e, e_end;
        tie(e, e_end) = in_edges(v, g);
        if (e == e_end)
            throw GraphException("Invalid hierarchical tree: No path from source to target.");
        v = source(*e, g);
        s_root.push_back(v);

        tie(e, e_end) = in_edges(u, g);
        if (e == e_end)
            throw GraphException("Invalid hierarchical tree: No path from source to target.");
        u = source(*e, g);
        if (u != v)
            t_root.push_back(u);
    }
    path = s_root;
165
    std::copy(t_root.rbegin(), t_root.rend(), std::back_inserter(path));
166
167
168
}


169
170
171
172
template <class Graph>
void graph_path(Graph& g, size_t s, size_t t, vector<size_t>& path)
{
    typename property_map_type::apply<size_t,
173
                                      typename property_map<Graph, vertex_index_t>::type>::type
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
        cpred;
    auto pred = cpred.get_unchecked(num_vertices(g));

    UndirectedAdaptor<Graph> ug(g);

    boost::breadth_first_search(ug, s,
                                boost::visitor(
                                    boost::make_bfs_visitor(
                                        boost::record_predecessors(
                                            pred,
                                            boost::on_tree_edge()))));
    size_t pos = t;
    path.push_back(pos);
    while (pos != s)
    {
        pos = pred[pos];
        path.push_back(pos);
    }
    std::reverse(path.begin(), path.end());
}

195
196
197
198
199
200
201
202
203
204
205
206
207
template<class T>
void pack(vector<point_t>& cp, vector<T>& ncp)
{
    ncp.resize(cp.size() * 2);
    for (size_t i = 0; i < cp.size(); ++i)
    {
        ncp[2 * i] = cp[i].first;
        ncp[2 * i + 1] = cp[i].second;
    }
}

struct do_get_cts
{
208
    template <class Graph, class Tree, class PosProp, class BProp, class CMap>
209
210
    void operator()(Graph& g, Tree* t, PosProp tpos, BProp beta, CMap cts,
                    bool is_tree, size_t max_depth) const
211
212
213
214
215
    {
        vector<size_t> path;
        vector<point_t> cp;
        vector<point_t> ncp;

216
        for (auto e : edges_range(g))
217
        {
218
219
            auto u = source(e, g);
            auto v = target(e, g);
220
221
            if (u == v)
                continue;
222

223
            path.clear();
224
            if (is_tree)
225
                tree_path(*t, u, v, path, max_depth);
226
227
            else
                graph_path(*t, u, v, path);
228
            cp.clear();
229
            get_control_points(path, tpos, beta[e], cp);
230
231
232
            ncp.clear();
            to_bezier(cp, ncp);
            transform(ncp);
233
            pack(ncp, cts[e]);
234
235
236
237
238
239
240
241
242
        }
    }
};

struct get_pointers
{
    template <class List>
    struct apply
    {
Tiago Peixoto's avatar
Tiago Peixoto committed
243
244
        typedef typename boost::mpl::transform<List,
                                               boost::mpl::quote1<std::add_pointer> >::type type;
245
246
247
    };
};

248
void get_cts(GraphInterface& gi, GraphInterface& tgi, boost::any otpos,
249
             boost::any obeta, boost::any octs, bool is_tree, size_t max_depth)
250
251
252
253
{
    typedef property_map_type::apply<vector<double>,
                                     GraphInterface::edge_index_map_t>::type
        eprop_t;
254
255
256
    typedef property_map_type::apply<double,
                                     GraphInterface::edge_index_map_t>::type
        beprop_t;
257
258

    eprop_t cts = boost::any_cast<eprop_t>(octs);
259
    beprop_t beta = boost::any_cast<beprop_t>(obeta);
260

261
    run_action<>()
Tiago Peixoto's avatar
Tiago Peixoto committed
262
        (gi, std::bind(do_get_cts(), placeholders::_1, placeholders::_2,
263
                       placeholders::_3, beta, cts, is_tree, max_depth),
264
265
         get_pointers::apply<graph_tool::detail::always_directed>::type(),
         vertex_scalar_vector_properties())
266
        (tgi.get_graph_view(), otpos);
267
}