root/library/bdm/base/bdmbase.cpp @ 624

Revision 624, 9.5 kB (checked in by smidl, 15 years ago)

RV accepts size=-1 as 'already existing size'

  • Property svn:eol-style set to native
Line 
1
2#include "bdmbase.h"
3
4//! Space of basic BDM structures
5namespace bdm {
6
7const int RV::BUFFER_STEP = 1;
8
9Array<string> RV::NAMES ( RV::BUFFER_STEP );
10
11ivec RV::SIZES ( RV::BUFFER_STEP );
12
13RV::str2int_map RV::MAP;
14
15void RV::clear_all() {
16        MAP.clear();
17        SIZES.clear();
18        NAMES = Array<string> ( BUFFER_STEP );
19}
20
21string RV::show_all(){
22  ostringstream os;
23  for(str2int_map::const_iterator iter=MAP.begin(); iter!=MAP.end(); iter++){
24          os << "key: " << iter->first << " val: " << iter->second <<endl;
25  }
26  return os.str();
27};
28
29int RV::init ( const string &name, int size ) {
30        //Refer
31        int id;
32        str2int_map::const_iterator iter = MAP.find ( name );
33        if ( iter == MAP.end() ) {
34                id = MAP.size() + 1;
35                //debug
36                /*              {
37                                        cout << endl;
38                                        str2int_map::const_iterator iter = MAP.begin();
39                                        for(str2int_map::const_iterator iter=MAP.begin(); iter!=MAP.end(); iter++){
40                                                cout << "key: " << iter->first << " val: " << iter->second <<endl;
41                                        }
42                                }*/
43
44                MAP.insert ( make_pair ( name, id ) ); //add new rv
45                if ( id >= NAMES.length() ) {
46                        NAMES.set_length ( id + BUFFER_STEP, true );
47                        SIZES.set_length ( id + BUFFER_STEP, true );
48                }
49                NAMES ( id ) = name;
50                SIZES ( id ) = size;
51                bdm_assert(size>0, "RV "+ name +" does not exists. Default size (-1) can not be assigned ");
52        } else {
53                id = iter->second;
54                if (size>0){
55                        bdm_assert ( SIZES ( id ) == size, "RV " + name + " of different size already exists" );
56                }
57        }
58        return id;
59};
60
61int RV::countsize() const {
62        int tmp = 0;
63        for ( int i = 0; i < len; i++ ) {
64                tmp += SIZES ( ids ( i ) );
65        }
66        return tmp;
67}
68
69ivec RV::cumsizes() const {
70        ivec szs ( len );
71        int tmp = 0;
72        for ( int i = 0; i < len; i++ ) {
73                tmp += SIZES ( ids ( i ) );
74                szs ( i ) = tmp;
75        }
76        return szs;
77}
78
79void RV::init ( const Array<std::string> &in_names, const ivec &in_sizes, const ivec &in_times ) {
80        len = in_names.length();
81        bdm_assert ( in_names.length() == in_times.length(), "check \"times\" " );
82        bdm_assert ( in_names.length() == in_sizes.length(), "check \"sizes\" " );
83
84        times.set_length ( len );
85        ids.set_length ( len );
86        int id;
87        for ( int i = 0; i < len; i++ ) {
88                id = init ( in_names ( i ), in_sizes ( i ) );
89                ids ( i ) = id;
90        }
91        times = in_times;
92        dsize = countsize();
93}
94
95RV::RV ( string name, int sz, int tm ) {
96        Array<string> A ( 1 );
97        A ( 0 ) = name;
98        init ( A, vec_1 ( sz ), vec_1 ( tm ) );
99}
100
101bool RV::add ( const RV &rv2 ) {
102        // TODO
103        if ( rv2.len > 0 ) { //rv2 is nonempty
104                ivec ind = rv2.findself ( *this ); //should be -1 all the time
105                ivec index = itpp::find ( ind == -1 );
106
107                if ( index.length() < rv2.len ) { //conflict
108                        ids = concat ( ids, rv2.ids ( index ) );
109                        times = concat ( times, rv2.times ( index ) );
110                } else {
111                        ids = concat ( ids, rv2.ids );
112                        times = concat ( times, rv2.times );
113                }
114                len = ids.length();
115                dsize = countsize();
116                return ( index.length() == rv2.len ); //conflict or not
117        } else { //rv2 is empty
118                return true; // no conflict
119        }
120};
121
122RV RV::subselect ( const ivec &ind ) const {
123        RV ret;
124        ret.ids = ids ( ind );
125        ret.times = times ( ind );
126        ret.len = ind.length();
127        ret.dsize = ret.countsize();
128        return ret;
129}
130
131RV RV::operator() ( int di1, int di2 ) const {
132        ivec sz = cumsizes();
133        int i1 = 0;
134        while ( sz ( i1 ) < di1 ) i1++;
135        int i2 = i1;
136        while ( sz ( i2 ) < di2 ) i2++;
137        return subselect ( linspace ( i1, i2 ) );
138}
139
140void RV::t_plus ( int delta ) {
141        times += delta;
142}
143
144bool RV::equal ( const RV &rv2 ) const {
145        return ( ids == rv2.ids ) && ( times == rv2.times );
146}
147
148shared_ptr<mpdf> epdf::condition ( const RV &rv ) const {
149        bdm_warning ( "Not implemented" );
150        return shared_ptr<mpdf>();
151}
152
153shared_ptr<epdf> epdf::marginal ( const RV &rv ) const {
154        bdm_warning ( "Not implemented" );
155        return shared_ptr<epdf>();
156}
157
158mat epdf::sample_m ( int N ) const {
159        mat X = zeros ( dim, N );
160        for ( int i = 0; i < N; i++ ) X.set_col ( i, this->sample() );
161        return X;
162}
163
164vec epdf::evallog_m ( const mat &Val ) const {
165        vec x ( Val.cols() );
166        for ( int i = 0; i < Val.cols(); i++ ) {
167                x ( i ) = evallog ( Val.get_col ( i ) );
168        }
169
170        return x;
171}
172
173vec epdf::evallog_m ( const Array<vec> &Avec ) const {
174        vec x ( Avec.size() );
175        for ( int i = 0; i < Avec.size(); i++ ) {
176                x ( i ) = evallog ( Avec ( i ) );
177        }
178
179        return x;
180}
181
182mat mpdf::samplecond_m ( const vec &cond, int N ) {
183        mat M ( dimension(), N );
184        for ( int i = 0; i < N; i++ ) {
185                M.set_col ( i, samplecond ( cond ) );
186        }
187
188        return M;
189}
190
191void mpdf::from_setting ( const Setting &set ) {
192        shared_ptr<RV> r = UI::build<RV> ( set, "rv", UI::optional );
193        if ( r ) {
194                set_rv ( *r );
195        }
196
197        r = UI::build<RV> ( set, "rvc", UI::optional );
198        if ( r ) {
199                set_rvc ( *r );
200        }
201}
202
203void datalink::set_connection ( const RV &rv, const RV &rv_up ) {
204        downsize = rv._dsize();
205        upsize = rv_up._dsize();
206        v2v_up = rv.dataind ( rv_up );
207        bdm_assert_debug ( v2v_up.length() == downsize, "rv is not fully in rv_up" );
208}
209
210void datalink::set_connection ( int ds, int us, const ivec &upind ) {
211        downsize = ds;
212        upsize = us;
213        v2v_up = upind;
214        bdm_assert_debug ( v2v_up.length() == downsize, "rv is not fully in rv_up" );
215}
216
217void datalink_part::set_connection ( const RV &rv, const RV &rv_up ) {
218        rv.dataind ( rv_up, v2v_down, v2v_up );
219        downsize = v2v_down.length();
220        upsize = v2v_up.length();
221}
222
223void datalink_m2e::set_connection ( const RV &rv, const RV &rvc, const RV &rv_up ) {
224        datalink::set_connection ( rv, rv_up );
225        condsize = rvc._dsize();
226        //establish v2c connection
227        rvc.dataind ( rv_up, v2c_lo, v2c_up );
228}
229
230vec datalink_m2e::get_cond ( const vec &val_up ) {
231        vec tmp ( condsize );
232        set_subvector ( tmp, v2c_lo, val_up ( v2c_up ) );
233        return tmp;
234}
235
236void datalink_m2e::pushup_cond ( vec &val_up, const vec &val, const vec &cond ) {
237        bdm_assert_debug ( downsize == val.length(), "Wrong val" );
238        bdm_assert_debug ( upsize == val_up.length(), "Wrong val_up" );
239        set_subvector ( val_up, v2v_up, val );
240        set_subvector ( val_up, v2c_up, cond );
241}
242
243std::ostream &operator<< ( std::ostream &os, const RV &rv ) {
244        int id;
245        for ( int i = 0; i < rv.len ; i++ ) {
246                id = rv.ids ( i );
247                os << id << "(" << RV::SIZES ( id ) << ")" <<  // id(size)=
248                "=" << RV::NAMES ( id )  << "_{"  << rv.times ( i ) << "}; "; //name_{time}
249        }
250        return os;
251}
252
253str RV::tostr() const {
254        ivec idlist ( dsize );
255        ivec tmlist ( dsize );
256        int i;
257        int pos = 0;
258        for ( i = 0; i < len; i++ ) {
259                idlist.set_subvector ( pos, pos + size ( i ) - 1, ids ( i ) );
260                tmlist.set_subvector ( pos, pos + size ( i ) - 1, times ( i ) );
261                pos += size ( i );
262        }
263        return str ( idlist, tmlist );
264}
265
266ivec RV::dataind ( const RV &rv2 ) const {
267        ivec res ( 0 );
268        if ( rv2._dsize() > 0 ) {
269                str str2 = rv2.tostr();
270                ivec part;
271                int i;
272                for ( i = 0; i < len; i++ ) {
273                        part = itpp::find ( ( str2.ids == ids ( i ) ) & ( str2.times == times ( i ) ) );
274                        res = concat ( res, part );
275                }
276        }
277
278        bdm_assert_debug ( res.length() == dsize, "this rv is not fully present in crv!" );
279        return res;
280
281}
282
283void RV::dataind ( const RV &rv2, ivec &selfi, ivec &rv2i ) const {
284        //clean results
285        selfi.set_size ( 0 );
286        rv2i.set_size ( 0 );
287
288        // just in case any rv is empty
289        if ( ( len == 0 ) || ( rv2.length() == 0 ) ) {
290                return;
291        }
292
293        //find comon rv
294        ivec cids = itpp::find ( this->findself ( rv2 ) >= 0 );
295
296        // index of
297        if ( cids.length() > 0 ) {
298                str str1 = tostr();
299                str str2 = rv2.tostr();
300
301                ivec part1;
302                ivec part2;
303                int i, j;
304                // find common rv in strs
305                for ( j = 0; j < cids.length(); j++ ) {
306                        i = cids ( j );
307                        part1 = itpp::find ( ( str1.ids == ids ( i ) ) & ( str1.times == times ( i ) ) );
308                        part2 = itpp::find ( ( str2.ids == ids ( i ) ) & ( str2.times == times ( i ) ) );
309                        selfi = concat ( selfi, part1 );
310                        rv2i = concat ( rv2i, part2 );
311                }
312        }
313        bdm_assert_debug ( selfi.length() == rv2i.length(), "this should not happen!" );
314}
315
316RV RV::subt ( const RV &rv2 ) const {
317        ivec res = this->findself ( rv2 ); // nonzeros
318        ivec valid;
319        if ( dsize > 0 ) {
320                valid = itpp::find ( res == -1 );    //-1 => value not found => it remains
321        }
322        return ( *this ) ( valid ); //keep those that were not found in rv2
323}
324
325ivec RV::findself ( const RV &rv2 ) const {
326        int i, j;
327        ivec tmp = -ones_i ( len );
328        for ( i = 0; i < len; i++ ) {
329                for ( j = 0; j < rv2.length(); j++ ) {
330                        if ( ( ids ( i ) == rv2.ids ( j ) ) & ( times ( i ) == rv2.times ( j ) ) ) {
331                                tmp ( i ) = j;
332                                break;
333                        }
334                }
335        }
336        return tmp;
337}
338
339ivec RV::findself_ids ( const RV &rv2 ) const {
340        int i, j;
341        ivec tmp = -ones_i ( len );
342        for ( i = 0; i < len; i++ ) {
343                for ( j = 0; j < rv2.length(); j++ ) {
344                        if ( ( ids ( i ) == rv2.ids ( j ) ) ) {
345                                tmp ( i ) = j;
346                                break;
347                        }
348                }
349        }
350        return tmp;
351}
352
353void RV::from_setting ( const Setting &set ) {
354        Array<string> A;
355        UI::get ( A, set, "names" );
356
357        ivec szs;
358        if ( !UI::get ( szs, set, "sizes" ) )
359                szs = ones_i ( A.length() );
360
361        ivec tms;
362        if ( !UI::get ( tms, set, "times" ) )
363                tms = zeros_i ( A.length() );
364
365        // TODO tady se bude plnit primo do jeho promennych, a pak se zavola validacnni metoda, takze cele prepsat, ano?
366        init ( A, szs, tms );
367}
368
369RV concat ( const RV &rv1, const RV &rv2 ) {
370        RV pom = rv1;
371        pom.add ( rv2 );
372        return pom;
373}
374
375void mepdf::from_setting ( const Setting &set ) {
376        shared_ptr<epdf> e ( UI::build<epdf> ( set, "epdf", UI::compulsory ) );
377        iepdf = e;
378        set_ep ( iepdf.get() );
379}
380
381RV get_composite_rv ( const Array<shared_ptr<mpdf> > &mpdfs,
382                      bool checkoverlap ) {
383        RV rv; //empty rv
384        bool rvaddok;
385        for ( int i = 0; i < mpdfs.length(); i++ ) {
386                rvaddok = rv.add ( mpdfs ( i )->_rv() ); //add rv to common rvs.
387                // If rvaddok==false, mpdfs overlap => assert error.
388                bdm_assert_debug ( rvaddok || !checkoverlap, "mprod::mprod() input mpdfs overlap in rv!" );
389        }
390
391        return rv;
392}
393
394void BM::bayesB ( const mat &Data ) {
395        for ( int t = 0; t < Data.cols(); t++ ) {
396                bayes ( Data.get_col ( t ) );
397        }
398}
399}
Note: See TracBrowser for help on using the browser.