source: issm/oecreview/Archive/12281-12300/ISSM-12291-12292.diff

Last change on this file was 12325, checked in by Eric.Larour, 13 years ago

11990 to 12321 oec compliance
File size: 15.2 KB
RevLine 
[12325]1Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/Container/Options.h
2===================================================================
3--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/Container/Options.h (revision 12291)
4+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/Container/Options.h (revision 12292)
5@@ -21,6 +21,8 @@
6 Option* GetOption(const char* name);
7 void Get(double* pvalue,const char* name);
8 void Get(double* pvalue,const char* name,double default_value);
9+ void Get(int* pvalue,const char* name);
10+ void Get(int* pvalue,const char* name,int default_value);
11 void Get(bool* pvalue,const char* name);
12 void Get(bool* pvalue,const char* name,bool default_value);
13 void Get(char** pvalue,const char* name);
14Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/Container/Observations.cpp
15===================================================================
16--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/Container/Observations.cpp (revision 12291)
17+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/Container/Observations.cpp (revision 12292)
18@@ -93,49 +93,73 @@
19
20 /*Methods*/
21 /*FUNCTION Observations::ObservationList{{{*/
22-void Observations::ObservationList(double **px,double **py,double **pobs,int* pnobs,double x_interp,double y_interp,double range){
23+void Observations::ObservationList(double **px,double **py,double **pobs,int* pnobs,double x_interp,double y_interp,double radius,int maxdata){
24
25 /*Output and Intermediaries*/
26- int nobs,i,index;
27+ bool stop;
28+ int nobs,i,j,k,n,counter;
29+ double h2,radius2;
30+ int *indices = NULL;
31+ double *dists = NULL;
32 double *x = NULL;
33 double *y = NULL;
34 double *obs = NULL;
35 Observation *observation = NULL;
36- int *indices = NULL;
37
38- /*Treat range*/
39- if(range==0){
40- range=this->quadtree->root->length;
41- }
42+ /*Compute radius square*/
43+ if(radius==0) radius=this->quadtree->root->length;
44+ radius2 = radius*radius;
45
46- /*Find all observations that are in range*/
47- nobs=0;
48- while(nobs==0){
49- xfree((void**)&indices);
50- this->quadtree->RangeSearch(&indices,&nobs,x_interp,y_interp,range);
51- if(!nobs){
52- /*No observation found, double range*/
53- range=2*range;
54+ /*Find all observations that are in radius*/
55+ indices = (int*)xmalloc(this->Size()*sizeof(int));
56+ dists = (double*)xmalloc(this->Size()*sizeof(double));
57+ nobs = 0;
58+
59+ for (i=0;i<this->Size();i++){
60+ observation=(Observation*)this->GetObjectByOffset(i);
61+ h2 = (observation->x-x_interp)*(observation->x-x_interp) + (observation->y-y_interp)*(observation->y-y_interp);
62+
63+ if(nobs==maxdata && h2>radius2) continue;
64+ if(nobs<=maxdata){
65+ indices[nobs] = i;
66+ dists[nobs] = h2;
67+ nobs++;
68 }
69- }
70- /*Check that we do not have more than 50 observations*/
71- //while(nobs>50){
72- // range=range/2;
73- // xfree((void**)&indices);
74- // this->quadtree->RangeSearch(&indices,&nobs,x_interp,y_interp,range);
75- //}
76+ if(nobs==1) continue;
77
78- if(!nobs) _error_("Unexpected error: no observation within current range and more than 50 observations if doubling range");
79+ /*Sort all dists up to now*/
80+ n=nobs-1;
81+ stop = false;
82+ for(k=0;k<n-1;k++){
83+ if(h2<dists[k]){
84+ counter=1;
85+ for(int jj=k;jj<n;jj++){
86+ j = n-counter;
87+ dists[j+1] = dists[j];
88+ indices[j+1] = indices[j];
89+ counter++;
90+ }
91+ dists[k] = h2;
92+ indices[k] = i;
93+ stop = true;
94+ break;
95+ }
96+ if(stop) break;
97+ }
98+ }
99+ xfree((void**)&dists);
100
101- /*Allocate vectors*/
102- x = (double*)xmalloc(nobs*sizeof(double));
103- y = (double*)xmalloc(nobs*sizeof(double));
104- obs = (double*)xmalloc(nobs*sizeof(double));
105+ if(nobs){
106+ /*Allocate vectors*/
107+ x = (double*)xmalloc(nobs*sizeof(double));
108+ y = (double*)xmalloc(nobs*sizeof(double));
109+ obs = (double*)xmalloc(nobs*sizeof(double));
110
111- /*Loop over all observations and fill in x, y and obs*/
112- for (i=0;i<nobs;i++){
113- observation=(Observation*)this->GetObjectByOffset(indices[i]);
114- observation->WriteXYObs(&x[i],&y[i],&obs[i]);
115+ /*Loop over all observations and fill in x, y and obs*/
116+ for (i=0;i<nobs;i++){
117+ observation=(Observation*)this->GetObjectByOffset(indices[i]);
118+ observation->WriteXYObs(&x[i],&y[i],&obs[i]);
119+ }
120 }
121
122 /*Assign output pointer*/
123Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/Container/Observations.h
124===================================================================
125--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/Container/Observations.h (revision 12291)
126+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/Container/Observations.h (revision 12292)
127@@ -22,7 +22,7 @@
128 ~Observations();
129
130 /*Methods*/
131- void ObservationList(double **px,double **py,double **pobs,int* pnobs,double x_interp,double y_interp,double range);
132+ void ObservationList(double **px,double **py,double **pobs,int* pnobs,double x_interp,double y_interp,double radius,int maxdata);
133 void QuadtreeColoring(double* A,double *x,double *y,int n);
134 void Variomap(double* gamma,double *x,int n);
135
136Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/Container/Options.cpp
137===================================================================
138--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/Container/Options.cpp (revision 12291)
139+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/Container/Options.cpp (revision 12292)
140@@ -66,6 +66,44 @@
141 return 1;
142 }
143 /*}}}*/
144+/*FUNCTION Options::Get(int* pvalue, char* name){{{1*/
145+void Options::Get(int* pvalue,const char* name){
146+
147+ vector<Object*>::iterator object;
148+ Option* option=NULL;
149+
150+ /*Get option*/
151+ option=GetOption(name);
152+
153+ /*If the pointer is not NULL, the option has been found*/
154+ if(option){
155+ option->Get(pvalue);
156+ }
157+ /*Else, the Option does not exist, no default provided*/
158+ else{
159+ _error_("option of name \"%s\" not found, and no default value has been provided",name);
160+ }
161+}
162+/*}}}*/
163+/*FUNCTION Options::Get(int* pvalue, char* name,int default_value){{{1*/
164+void Options::Get(int* pvalue,const char* name,int default_value){
165+
166+ vector<Object*>::iterator object;
167+ Option* option=NULL;
168+
169+ /*Get option*/
170+ option=GetOption(name);
171+
172+ /*If the pointer is not NULL, the option has been found*/
173+ if(option){
174+ option->Get(pvalue);
175+ }
176+ /*Else, the Option does not exist, a default is provided here*/
177+ else{
178+ *pvalue=default_value;
179+ }
180+}
181+/*}}}*/
182 /*FUNCTION Options::Get(double* pvalue, char* name){{{1*/
183 void Options::Get(double* pvalue,const char* name){
184
185Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/shared/Threads/LaunchThread.cpp
186===================================================================
187--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/shared/Threads/LaunchThread.cpp (revision 12291)
188+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/shared/Threads/LaunchThread.cpp (revision 12292)
189@@ -66,5 +66,4 @@
190
191 function((void*)&handle);
192 #endif
193-
194 }
195Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/modules/Krigingx/Krigingx.cpp
196===================================================================
197--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/modules/Krigingx/Krigingx.cpp (revision 12291)
198+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/modules/Krigingx/Krigingx.cpp (revision 12292)
199@@ -23,7 +23,8 @@
200 double *error = NULL;
201
202 /*Intermediaries*/
203- double range;
204+ int mindata,maxdata;
205+ double radius;
206 char *output = NULL;
207 Variogram *variogram = NULL;
208 Observations *observations = NULL;
209@@ -37,7 +38,9 @@
210
211 /*Get Variogram from Options*/
212 ProcessVariogram(&variogram,options);
213- options->Get(&range,"searchrange",0.);
214+ options->Get(&radius,"searchradius",0.);
215+ options->Get(&mindata,"mindata",1);
216+ options->Get(&maxdata,"maxdata",50);
217
218 /*Process observation dataset*/
219 observations=new Observations(obs_list,obs_x,obs_y,obs_length,options);
220@@ -61,7 +64,9 @@
221 gate.n_interp = n_interp;
222 gate.x_interp = x_interp;
223 gate.y_interp = y_interp;
224- gate.range = range;
225+ gate.radius = radius;
226+ gate.mindata = mindata;
227+ gate.maxdata = maxdata;
228 gate.variogram = variogram;
229 gate.observations = observations;
230 gate.predictions = predictions;
231@@ -105,7 +110,9 @@
232 int n_interp = gate->n_interp;
233 double *x_interp = gate->x_interp;
234 double *y_interp = gate->y_interp;
235- double range = gate->range;
236+ double radius = gate->radius;
237+ int mindata = gate->mindata;
238+ int maxdata = gate->maxdata;
239 Variogram *variogram = gate->variogram;
240 Observations *observations = gate->observations;
241 double *predictions = gate->predictions;
242@@ -136,10 +143,10 @@
243 printf("\r interpolation progress: %5.2lf%%",localpercent);
244
245 /*Get list of observations for current point*/
246- observations->ObservationList(&x,&y,&obs,&n_obs,x_interp[idx],y_interp[idx],range);
247- if(n_obs==0){
248- predictions[idx] = -999.;
249- error[idx] = -999.;
250+ observations->ObservationList(&x,&y,&obs,&n_obs,x_interp[idx],y_interp[idx],radius,maxdata);
251+ if(n_obs<mindata){
252+ predictions[idx] = -999.0;
253+ error[idx] = -999.0;
254 continue;
255 }
256
257Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/modules/Krigingx/Krigingx.h
258===================================================================
259--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/modules/Krigingx/Krigingx.h (revision 12291)
260+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/modules/Krigingx/Krigingx.h (revision 12292)
261@@ -20,7 +20,9 @@
262 int n_interp;
263 double *x_interp;
264 double *y_interp;
265- double range;
266+ double radius;
267+ int mindata;
268+ int maxdata;
269 Variogram *variogram;
270 Observations *observations;
271 double *predictions;
272Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/OptionStruct.h
273===================================================================
274--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/OptionStruct.h (revision 12291)
275+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/OptionStruct.h (revision 12292)
276@@ -38,6 +38,7 @@
277 int NumEl();
278 int NDims();
279 int* Size();
280+ void Get(int* pvalue){_error_("An OptionStruct object cannot return a int");};
281 void Get(double* pvalue){_error_("An OptionStruct object cannot return a double");};
282 void Get(bool* pvalue){ _error_("An OptionStruct object cannot return a bool");};
283 void Get(char** pvalue){ _error_("An OptionStruct object cannot return a string");};
284Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/OptionDouble.cpp
285===================================================================
286--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/OptionDouble.cpp (revision 12291)
287+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/OptionDouble.cpp (revision 12292)
288@@ -119,6 +119,18 @@
289 return(Option::Size());
290 }
291 /*}}}*/
292+/*FUNCTION OptionDouble::Get(int* pvalue) {{{1*/
293+void OptionDouble::Get(int* pvalue){
294+
295+ /*We should first check that the size is one*/
296+ if(this->NumEl()!=1){
297+ _error_("option \"%s\" has %i elements and cannot return a single int",this->name,this->NumEl());
298+ }
299+
300+ /*Assign output pointer*/
301+ *pvalue=(int)this->values[0];
302+}
303+/*}}}*/
304 /*FUNCTION OptionDouble::Get(double* pvalue) {{{1*/
305 void OptionDouble::Get(double* pvalue){
306
307Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/OptionDouble.h
308===================================================================
309--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/OptionDouble.h (revision 12291)
310+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/OptionDouble.h (revision 12292)
311@@ -38,6 +38,7 @@
312 int NumEl();
313 int NDims();
314 int* Size();
315+ void Get(int* pvalue);
316 void Get(double* pvalue);
317 void Get(bool* pvalue){ _error_("An OptionDouble object cannot return a bool");};
318 void Get(char** pvalue){ _error_("An OptionDouble object cannot return a string");};
319Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/OptionLogical.h
320===================================================================
321--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/OptionLogical.h (revision 12291)
322+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/OptionLogical.h (revision 12292)
323@@ -38,6 +38,7 @@
324 int NumEl();
325 int NDims();
326 int* Size();
327+ void Get(int* pvalue){_error_("An OptionLogical object cannot return a int");};
328 void Get(double* pvalue){_error_("An OptionLogical object cannot return a double");};
329 void Get(bool* pvalue);
330 void Get(char** pvalue){ _error_("An OptionLogical object cannot return a string");};
331Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/OptionChar.h
332===================================================================
333--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/OptionChar.h (revision 12291)
334+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/OptionChar.h (revision 12292)
335@@ -38,6 +38,7 @@
336 int NumEl();
337 int NDims();
338 int* Size();
339+ void Get(int* pvalue){_error_("An OptionChar object cannot return a int");};
340 void Get(double* pvalue){_error_("An OptionChar object cannot return a double");};
341 void Get(bool* pvalue){ _error_("An OptionChar object cannot return a bool");};
342 void Get(char** pvalue);
343Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/Option.h
344===================================================================
345--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/Option.h (revision 12291)
346+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/Option.h (revision 12292)
347@@ -41,6 +41,7 @@
348 virtual int NumEl()=0;
349 virtual int NDims()=0;
350 virtual int* Size()=0;
351+ virtual void Get(int* pvalue)=0;
352 virtual void Get(double* pvalue)=0;
353 virtual void Get(bool* pvalue)=0;
354 virtual void Get(char** pvalue)=0;
355Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/OptionCell.h
356===================================================================
357--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/OptionCell.h (revision 12291)
358+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Options/OptionCell.h (revision 12292)
359@@ -38,6 +38,7 @@
360 int NumEl();
361 int NDims();
362 int* Size();
363+ void Get(int* pvalue){_error_("An OptionCell object cannot return a int");};
364 void Get(double* pvalue){_error_("An OptionCell object cannot return a double");};
365 void Get(bool* pvalue){ _error_("An OptionCell object cannot return a bool");};
366 void Get(char** pvalue){ _error_("An OptionCell object cannot return a string");};
367Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Kriging/GaussianVariogram.cpp
368===================================================================
369--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Kriging/GaussianVariogram.cpp (revision 12291)
370+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Kriging/GaussianVariogram.cpp (revision 12292)
371@@ -74,7 +74,7 @@
372 a = 1./3.;
373 cova = (sill-nugget)*exp(-h2/(a*range*range));
374
375- //if(h2<0.0000001) cova = sill;
376+ if(h2<0.0000001) cova = sill;
377
378 return cova;
379 }
Note: See TracBrowser for help on using the repository browser.