Changeset 1395 for applications/robust/robustlib.cpp

Timestamp:

09/21/11 20:35:02 (13 years ago)

Author:

sindj

Message:

Dodelavani toho zatracenyho samplingu. JS

Files:

• applications/robust/robustlib.cpp (modified) (4 diffs)

applications/robust/robustlib.cpp

@@ -210 +210 @@
                                         }
                                         
-                                        double factor_multiplier = pow((double)-1,k)/pow((*a_ref).first,sigma_order-k);//pow((double)-1,k+1)
+                                        double factor_multiplier = pow((double)-1,k)/pow((*a_ref).first,sigma_order-a_order+k+1);//pow((double)-1,k+1)
                                                                                 
+                                        double value = 0;
                                         if(k!=0)
-                                        {       
-                                                double value = 0;
+                                        {                                                       
                                                 ivec control_vec = ivec();
                                                 control_vec.ins(0,my_emlig->number_of_parameters-a_order+1);                                            
-                                                for(multiset<my_ivec>::iterator combi_ref = this->my_emlig->correction_factors[k-1].begin();combi_ref!=this->my_emlig->correction_factors[k-1].upper_bound(my_ivec(control_vec));combi_ref++)
+                                                for(multiset<my_ivec>::iterator combi_ref = this->my_emlig->correction_factors[a_order-k-1].begin();combi_ref!=this->my_emlig->correction_factors[a_order-k-1].upper_bound(my_ivec(control_vec));combi_ref++)
                                                 {
-                                                        double bracket_factor = 1/fact(a_order-1-k);
+                                                        double bracket_factor = 1/fact(k);
                                                          
                                                         for(int j = 0;j<(*combi_ref).size();j++)
@@ -229 +229 @@
                                                                                                         
                                                         value += bracket_factor*factor_multiplier*k_multiplier;                                                 
-                                                }
-
-                                                simplex->insert_gamma(k,value,(*a_ref).first);
-                                                gamma_facs[k] += value;
+                                                }                                               
                                         }
                                         else
                                         {
-                                                double value = k_multiplier*factor_multiplier/fact(a_order-1);                                          
-                                                
-                                                simplex->insert_gamma(0,value,(*a_ref).first);
-                                                gamma_facs[0] += value;
-                                        }                                               
+                                                value = k_multiplier*factor_multiplier;                                                 
+                                        }
+
+                                        simplex->insert_gamma(k,value,(*a_ref).first);
+                                        gamma_facs[k] += value;
                                 }                               
                         }
 
                         double int_value = 0;
+                        double volume = det_jacobian/fact(jacobian.rows());
                         int gamma_size = (int)gamma_facs.size();
                         if(sigma_order-gamma_size>=0)
@@ -250 +248 @@
                                 for(int k = 0;k<gamma_size;k++)
                                 {       
-                                        int_value += det_jacobian*fact(sigma_order-gamma_size+k)/fact(jacobian.rows())*gamma_facs[k]/pow(2.0,sigma_order);                              
+                                        int_value += volume*fact(sigma_order+k)*gamma_facs[k]/pow(2.0,((toprow*)this)->condition_order);                                
                                 }       
                         }
@@ -259 +257 @@
                         
                         simplex->probability = int_value;
+                        simplex->volume = volume;
                         //cout << "Probability:" << int_value << endl;
+                        //pause(0.1);
                         return int_value;               
                 }