Changeset 1389 for applications/robust/robustlib.h

Timestamp:

09/13/11 19:44:06 (13 years ago)

Author:

sindj

Message:

Odstraneni chyby v samplingu. JS

Files:

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

applications/robust/robustlib.h

@@ -21 +21 @@
 using namespace std;
 using namespace itpp;
+
+static Exponential_RNG ExpRNG;
 
 const double max_range = 5;//numeric_limits<double>::max()/10e-10;
@@ -1372 +1374 @@
                 {
                         condition_to_add = new condition(toadd);                        
-                        conditions.push_back(new_condition);                    
+                        conditions.push_back(condition_to_add);                 
                 }               
                 
@@ -2618 +2620 @@
                         pair<double,double> probability_and_sigma = choose_sigma(condition_and_simplex.second);
 
-                        int dimension = condition_and_simplex.second->vertices.size()-1;
-
-                        mat jacobian(dimension,dimension);
-                        vec gradient = condition_and_simplex.first.right(dimension);
-
-                        vertex* base_vert = *condition_and_simplex.second->vertices.begin();
-
-                        //// cout << "Base vertex coords(should be close to est. param.): " << base_vert->get_coordinates() << endl;
+                        int dimension = condition_and_simplex.second->vertices.size();
+
+                        mat jacobian(dimension,dimension-1);
+                        vec gradient = condition_and_simplex.first;
                                 
                         int row_count = 0;
 
-                        for(set<vertex*>::iterator vert_ref = ++condition_and_simplex.second->vertices.begin();vert_ref!=condition_and_simplex.second->vertices.end();vert_ref++)
-                        {
-                                vec current_coords = (*vert_ref)->get_coordinates();
-
-                                //// cout << "Coords of vertex[" << row_count << "]: " << current_coords << endl; 
-                                        
-                                vec relative_coords = current_coords-base_vert->get_coordinates();                              
-
-                                jacobian.set_row(row_count,relative_coords);
-
+                        for(set<vertex*>::iterator vert_ref = condition_and_simplex.second->vertices.begin();vert_ref!=condition_and_simplex.second->vertices.end();vert_ref++)
+                        {
+                                jacobian.set_row(row_count,(*vert_ref)->get_coordinates());
                                 row_count++;
-                        }                               
-                                
-                        //// cout << "Jacobian: " << jacobian << endl;                  
-
-                        /// \todo Is this correct? Are the random coordinates really jointly uniform? I don't know.
+                        }               
+                        
+                        ExpRNG.setup(1);
+
                         vec sample_coords;
-                        double sampling_diff = 1;
-                        for(int j = 0;j<number_of_parameters;j++)
-                        {
-                                double rnumber = randu()*sampling_diff;
-
-                                sample_coords.ins(0,rnumber);
-
-                                sampling_diff -= rnumber;
+                        double sample_sum = 0;
+                        for(int j = 0;j<dimension;j++)
+                        {
+                                double rnumber = ExpRNG();
+
+                                sample_sum += rnumber;
+
+                                sample_coords.ins(0,rnumber);                           
                         }
 
-                        sample_coords = jacobian.T()*sample_coords+(*base_vert).get_coordinates();
+                        sample_coords /= sample_sum;
+
+                        sample_coords = jacobian.T()*sample_coords;
 
                         vec extended_coords = sample_coords;
@@ -2662 +2654 @@
 
                         double exponent = extended_coords*condition_and_simplex.first;
-                        double sample_prob = 1/condition_and_simplex.second->probability/pow(probability_and_sigma.second,(int)conditions.size()-number_of_parameters+3)*exp((-1)/probability_and_sigma.second*exponent);
+                        double sample_prob = 1/pow(probability_and_sigma.second,(int)conditions.size()-number_of_parameters+3)*exp((-1)/probability_and_sigma.second*exponent);
                         sample_prob *= probability_and_sigma.first;