Show
Ignore:
Timestamp:
08/18/08 14:27:50 (16 years ago)
Author:
smidl
Message:

Documentation update

Files:
1 modified

Legend:

Unmodified
Added
Removed
  • doc/html/libPF_8h-source.html

    r106 r139  
    55<link href="tabs.css" rel="stylesheet" type="text/css"> 
    66</head><body> 
    7 <!-- Generated by Doxygen 1.5.5 --> 
     7<!-- Generated by Doxygen 1.5.6 --> 
    88<div class="navigation" id="top"> 
    99  <div class="tabs"> 
     
    112112<a name="l00131"></a>00131  
    113113<a name="l00132"></a>00132 <span class="comment">//SimStr:</span> 
    114 <a name="l00133"></a>00133         <span class="keywordtype">double</span> SSAT;     
     114<a name="l00133"></a>00133         <span class="keywordtype">double</span> SSAT; 
    115115<a name="l00134"></a>00134 }; 
    116116<a name="l00135"></a>00135  
     
    128128<a name="l00147"></a>00147         <span class="keywordtype">double</span> WL = 0.0; 
    129129<a name="l00148"></a>00148  
    130 <a name="l00149"></a>00149         <span class="keywordflow">for</span> ( i=0;i&lt;<a class="code" href="classPF.html#2c2f44ed7a4eaa42e07bdb58d503f280" title="number of particles;">n</a>;i++ ) { 
    131 <a name="l00150"></a>00150                 <span class="comment">//generate new samples from paramater evolution model;</span> 
    132 <a name="l00151"></a>00151                 <a class="code" href="classPF.html#cf7dad75e31215780a746c30e71ad9c5" title="pointer into eEmp ">_samples</a> ( i ) = <a class="code" href="classPF.html#d92ac103f88f8c21e197e90af5695a09" title="Parameter evolution model.">par</a>.<a class="code" href="classmpdf.html#b0193a350c97933ddf15b15a130da352" title="Returns the required moment of the epdf.">samplecond</a> ( <a class="code" href="classPF.html#cf7dad75e31215780a746c30e71ad9c5" title="pointer into eEmp ">_samples</a> ( i ), llsP ( i ) ); 
    133 <a name="l00152"></a>00152                 Bms[i]-&gt;condition ( <a class="code" href="classPF.html#cf7dad75e31215780a746c30e71ad9c5" title="pointer into eEmp ">_samples</a> ( i ) ); 
    134 <a name="l00153"></a>00153                 Bms[i]-&gt;bayes ( dt ); 
    135 <a name="l00154"></a>00154                 lls ( i ) = Bms[i]-&gt;_ll(); <span class="comment">// lls above is also in proposal her must be lls(i) =, not +=!!</span> 
    136 <a name="l00155"></a>00155  
     130<a name="l00149"></a>00149 <span class="preprocessor">        #pragma omp parallel for</span> 
     131<a name="l00150"></a>00150 <span class="preprocessor"></span>        <span class="keywordflow">for</span> ( i=0;i&lt;<a class="code" href="classPF.html#2c2f44ed7a4eaa42e07bdb58d503f280" title="number of particles;">n</a>;i++ ) {    
     132<a name="l00151"></a>00151                 <span class="comment">//generate new samples from paramater evolution model;</span> 
     133<a name="l00152"></a>00152                 <a class="code" href="classPF.html#cf7dad75e31215780a746c30e71ad9c5" title="pointer into eEmp ">_samples</a> ( i ) = <a class="code" href="classPF.html#d92ac103f88f8c21e197e90af5695a09" title="Parameter evolution model.">par</a>.<a class="code" href="classmpdf.html#b0193a350c97933ddf15b15a130da352" title="Returns the required moment of the epdf.">samplecond</a> ( <a class="code" href="classPF.html#cf7dad75e31215780a746c30e71ad9c5" title="pointer into eEmp ">_samples</a> ( i ), llsP ( i ) ); 
     134<a name="l00153"></a>00153                 Bms[i]-&gt;condition ( <a class="code" href="classPF.html#cf7dad75e31215780a746c30e71ad9c5" title="pointer into eEmp ">_samples</a> ( i ) ); 
     135<a name="l00154"></a>00154                 Bms[i]-&gt;bayes ( dt ); 
     136<a name="l00155"></a>00155                 lls ( i ) = Bms[i]-&gt;_ll(); <span class="comment">// lls above is also in proposal her must be lls(i) =, not +=!!</span> 
    137137<a name="l00156"></a>00156                 <span class="keywordflow">if</span> ( lls ( i ) &gt;mlls ) mlls=lls ( i ); <span class="comment">//find maximum likelihood (for numerical stability)</span> 
    138138<a name="l00157"></a>00157         } 
    139 <a name="l00158"></a>00158          
    140 <a name="l00159"></a>00159         <span class="keywordflow">if</span> ( <span class="keyword">true</span> ) { 
    141 <a name="l00160"></a>00160                 <span class="keywordflow">for</span> ( i=0;i&lt;n;i++ ) { 
    142 <a name="l00161"></a>00161                         WL = <a class="code" href="classPF.html#5c87aba508df321ff26536ced64dbb3a" title="pointer into eEmp ">_w</a> ( i ) *exp ( llsP ( i ) ); <span class="comment">//using old weights!</span> 
    143 <a name="l00162"></a>00162                         sumLWL += exp ( lls ( i ) ) *WL; 
    144 <a name="l00163"></a>00163                         sumL2WL += exp ( 2*lls ( i ) ) *WL; 
    145 <a name="l00164"></a>00164                 } 
    146 <a name="l00165"></a>00165                 SSAT  = sumL2WL/(sumLWL*sumLWL); 
    147 <a name="l00166"></a>00166         } 
    148 <a name="l00167"></a>00167  
    149 <a name="l00168"></a>00168         <span class="comment">// compute weights</span> 
    150 <a name="l00169"></a>00169         <span class="keywordflow">for</span> ( i=0;i&lt;n;i++ ) { 
    151 <a name="l00170"></a>00170                 <a class="code" href="classPF.html#5c87aba508df321ff26536ced64dbb3a" title="pointer into eEmp ">_w</a> ( i ) *= exp ( lls ( i ) - mlls ); <span class="comment">// multiply w by likelihood</span> 
    152 <a name="l00171"></a>00171         } 
    153 <a name="l00172"></a>00172  
    154 <a name="l00173"></a>00173         <span class="keywordflow">if</span> (sum(<a class="code" href="classPF.html#5c87aba508df321ff26536ced64dbb3a" title="pointer into eEmp ">_w</a>)&gt;0.0){ 
    155 <a name="l00174"></a>00174                 <a class="code" href="classPF.html#5c87aba508df321ff26536ced64dbb3a" title="pointer into eEmp ">_w</a> /=sum ( <a class="code" href="classPF.html#5c87aba508df321ff26536ced64dbb3a" title="pointer into eEmp ">_w</a> ); <span class="comment">//?</span> 
    156 <a name="l00175"></a>00175                 } 
    157 <a name="l00176"></a>00176 <span class="keywordflow">else</span> 
    158 <a name="l00177"></a>00177 { 
    159 <a name="l00178"></a>00178 cout&lt;&lt;<span class="stringliteral">"sum(w)==0"</span>&lt;&lt;endl; 
    160 <a name="l00179"></a>00179 } 
    161 <a name="l00180"></a>00180  
    162 <a name="l00181"></a>00181  
    163 <a name="l00182"></a>00182         <span class="keywordtype">double</span> eff = 1.0/ ( <a class="code" href="classPF.html#5c87aba508df321ff26536ced64dbb3a" title="pointer into eEmp ">_w</a>*<a class="code" href="classPF.html#5c87aba508df321ff26536ced64dbb3a" title="pointer into eEmp ">_w</a> ); 
    164 <a name="l00183"></a>00183         <span class="keywordflow">if</span> ( eff &lt; ( 0.3*n ) ) { 
    165 <a name="l00184"></a>00184                 ind = <a class="code" href="classPF.html#1a0a09e309da997f63ae8e30d1e9806b" title="posterior density">est</a>.<a class="code" href="classeEmp.html#77268292fc4465cb73ddbfb1f2932a59" title="Function performs resampling, i.e. removal of low-weight samples and duplication...">resample</a>(); 
    166 <a name="l00185"></a>00185                 <span class="comment">// Resample Bms!</span> 
    167 <a name="l00186"></a>00186  
    168 <a name="l00187"></a>00187                 <span class="keywordflow">for</span> ( i=0;i&lt;n;i++ ) { 
    169 <a name="l00188"></a>00188                         <span class="keywordflow">if</span> ( ind ( i ) !=i ) {<span class="comment">//replace the current Bm by a new one</span> 
    170 <a name="l00189"></a>00189                                 <span class="comment">//fixme this would require new assignment operator</span> 
    171 <a name="l00190"></a>00190                                 <span class="comment">// *Bms[i] = *Bms[ind ( i ) ];</span> 
    172 <a name="l00191"></a>00191  
    173 <a name="l00192"></a>00192                                 <span class="comment">// poor-man's solution: replicate constructor here</span> 
    174 <a name="l00193"></a>00193                                 <span class="comment">// copied from MPF::MPF</span> 
    175 <a name="l00194"></a>00194                                 <span class="keyword">delete</span> Bms[i]; 
    176 <a name="l00195"></a>00195                                 Bms[i] = <span class="keyword">new</span> BM_T ( *Bms[ind ( i ) ] ); <span class="comment">//copy constructor</span> 
    177 <a name="l00196"></a>00196                                 <a class="code" href="classepdf.html" title="Probability density function with numerical statistics, e.g. posterior density.">epdf</a>&amp; pom=Bms[i]-&gt;_epdf(); 
    178 <a name="l00197"></a>00197                                 jest.set_elements ( i,1.0/n,&amp;pom ); 
    179 <a name="l00198"></a>00198                         } 
    180 <a name="l00199"></a>00199                 }; 
    181 <a name="l00200"></a>00200                 cout &lt;&lt; <span class="charliteral">'.'</span>; 
    182 <a name="l00201"></a>00201         } 
    183 <a name="l00202"></a>00202 } 
    184 <a name="l00203"></a>00203  
    185 <a name="l00204"></a>00204 <span class="preprocessor">#endif // KF_H</span> 
    186 <a name="l00205"></a>00205 <span class="preprocessor"></span> 
     139<a name="l00158"></a>00158  
     140<a name="l00159"></a>00159         <span class="keywordflow">if</span> ( <span class="keyword">false</span>) { 
     141<a name="l00160"></a>00160 <span class="preprocessor">                #pragma omp parallel for reduction(+:sumLWL,sumL2WL) private(WL)</span> 
     142<a name="l00161"></a>00161 <span class="preprocessor"></span>                <span class="keywordflow">for</span> ( i=0;i&lt;n;i++ ) { 
     143<a name="l00162"></a>00162                         WL = <a class="code" href="classPF.html#5c87aba508df321ff26536ced64dbb3a" title="pointer into eEmp ">_w</a> ( i ) *exp ( llsP ( i ) ); <span class="comment">//using old weights!</span> 
     144<a name="l00163"></a>00163                         sumLWL += exp ( lls ( i ) ) *WL; 
     145<a name="l00164"></a>00164                         sumL2WL += exp ( 2*lls ( i ) ) *WL; 
     146<a name="l00165"></a>00165                 } 
     147<a name="l00166"></a>00166                 SSAT  = sumL2WL/ ( sumLWL*sumLWL ); 
     148<a name="l00167"></a>00167         } 
     149<a name="l00168"></a>00168  
     150<a name="l00169"></a>00169         <span class="keywordtype">double</span> sum_w=0.0; 
     151<a name="l00170"></a>00170         <span class="comment">// compute weights</span> 
     152<a name="l00171"></a>00171 <span class="preprocessor">        #pragma omp parallel for</span> 
     153<a name="l00172"></a>00172 <span class="preprocessor"></span>        <span class="keywordflow">for</span> ( i=0;i&lt;n;i++ ) { 
     154<a name="l00173"></a>00173                 <a class="code" href="classPF.html#5c87aba508df321ff26536ced64dbb3a" title="pointer into eEmp ">_w</a> ( i ) *= exp ( lls ( i ) - mlls ); <span class="comment">// multiply w by likelihood</span> 
     155<a name="l00174"></a>00174                 sum_w+=<a class="code" href="classPF.html#5c87aba508df321ff26536ced64dbb3a" title="pointer into eEmp ">_w</a>(i); 
     156<a name="l00175"></a>00175         } 
     157<a name="l00176"></a>00176  
     158<a name="l00177"></a>00177         <span class="keywordflow">if</span> ( sum_w  &gt;0.0 ) { 
     159<a name="l00178"></a>00178                 <a class="code" href="classPF.html#5c87aba508df321ff26536ced64dbb3a" title="pointer into eEmp ">_w</a> /=sum_w; <span class="comment">//?</span> 
     160<a name="l00179"></a>00179         } <span class="keywordflow">else</span> { 
     161<a name="l00180"></a>00180                 cout&lt;&lt;<span class="stringliteral">"sum(w)==0"</span>&lt;&lt;endl; 
     162<a name="l00181"></a>00181         } 
     163<a name="l00182"></a>00182  
     164<a name="l00183"></a>00183  
     165<a name="l00184"></a>00184         <span class="keywordtype">double</span> eff = 1.0/ ( <a class="code" href="classPF.html#5c87aba508df321ff26536ced64dbb3a" title="pointer into eEmp ">_w</a>*<a class="code" href="classPF.html#5c87aba508df321ff26536ced64dbb3a" title="pointer into eEmp ">_w</a> ); 
     166<a name="l00185"></a>00185         <span class="keywordflow">if</span> ( eff &lt; ( 0.3*n ) ) { 
     167<a name="l00186"></a>00186                 ind = <a class="code" href="classPF.html#1a0a09e309da997f63ae8e30d1e9806b" title="posterior density">est</a>.<a class="code" href="classeEmp.html#77268292fc4465cb73ddbfb1f2932a59" title="Function performs resampling, i.e. removal of low-weight samples and duplication...">resample</a>(); 
     168<a name="l00187"></a>00187                 <span class="comment">// Resample Bms!</span> 
     169<a name="l00188"></a>00188  
     170<a name="l00189"></a>00189 <span class="preprocessor">                #pragma omp parallel for</span> 
     171<a name="l00190"></a>00190 <span class="preprocessor"></span>                <span class="keywordflow">for</span> ( i=0;i&lt;n;i++ ) { 
     172<a name="l00191"></a>00191                         <span class="keywordflow">if</span> ( ind ( i ) !=i ) {<span class="comment">//replace the current Bm by a new one</span> 
     173<a name="l00192"></a>00192                                 <span class="comment">//fixme this would require new assignment operator</span> 
     174<a name="l00193"></a>00193                                 <span class="comment">// *Bms[i] = *Bms[ind ( i ) ];</span> 
     175<a name="l00194"></a>00194  
     176<a name="l00195"></a>00195                                 <span class="comment">// poor-man's solution: replicate constructor here</span> 
     177<a name="l00196"></a>00196                                 <span class="comment">// copied from MPF::MPF</span> 
     178<a name="l00197"></a>00197                                 <span class="keyword">delete</span> Bms[i]; 
     179<a name="l00198"></a>00198                                 Bms[i] = <span class="keyword">new</span> BM_T ( *Bms[ind ( i ) ] ); <span class="comment">//copy constructor</span> 
     180<a name="l00199"></a>00199                                 <a class="code" href="classepdf.html" title="Probability density function with numerical statistics, e.g. posterior density.">epdf</a>&amp; pom=Bms[i]-&gt;_epdf(); 
     181<a name="l00200"></a>00200                                 jest.set_elements ( i,1.0/n,&amp;pom ); 
     182<a name="l00201"></a>00201                         } 
     183<a name="l00202"></a>00202                 }; 
     184<a name="l00203"></a>00203                 cout &lt;&lt; <span class="charliteral">'.'</span>; 
     185<a name="l00204"></a>00204         } 
     186<a name="l00205"></a>00205 } 
    187187<a name="l00206"></a>00206  
     188<a name="l00207"></a>00207 <span class="preprocessor">#endif // KF_H</span> 
     189<a name="l00208"></a>00208 <span class="preprocessor"></span> 
     190<a name="l00209"></a>00209  
    188191</pre></div></div> 
    189 <hr size="1"><address style="text-align: right;"><small>Generated on Fri May 9 23:06:27 2008 for mixpp by&nbsp; 
     192<hr size="1"><address style="text-align: right;"><small>Generated on Sat Aug 16 11:58:41 2008 for mixpp by&nbsp; 
    190193<a href="http://www.doxygen.org/index.html"> 
    191 <img src="doxygen.png" alt="doxygen" align="middle" border="0"></a> 1.5.5 </small></address> 
     194<img src="doxygen.png" alt="doxygen" align="middle" border="0"></a> 1.5.6 </small></address> 
    192195</body> 
    193196</html>