source: issm/oecreview/Archive/12221-12240/ISSM-12230-12231.diff@ 12325

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

11990 to 12321 oec compliance
File size: 8.3 KB
RevLine 
[12325]1Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/Container/Observations.cpp
2===================================================================
3--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/Container/Observations.cpp (revision 12230)
4+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/Container/Observations.cpp (revision 12231)
5@@ -112,12 +112,12 @@
6
7 if(nobs==0){
8 /*No observation found, double range*/
9- printf("No observation found within range, doubling range\n");
10+ //printf("No observation found within range, doubling range\n");
11 xfree((void**)&indices);
12 this->ObservationList(&x,&y,&obs,&nobs,x_interp,y_interp,range*2);
13 }
14 else{
15- if(nobs>1000) printf("Taking more than 1000 observations\n");
16+ //if(nobs>1000) printf("Taking more than 1000 observations\n");
17 /*Allocate vectors*/
18 x = (double*)xmalloc(nobs*sizeof(double));
19 y = (double*)xmalloc(nobs*sizeof(double));
20Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/modules/Krigingx/Krigingx.cpp
21===================================================================
22--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/modules/Krigingx/Krigingx.cpp (revision 12230)
23+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/modules/Krigingx/Krigingx.cpp (revision 12231)
24@@ -15,28 +15,25 @@
25 #endif
26
27 #include "../../objects/Kriging/GaussianVariogram.h"
28+/*FUNCTION Krigingx{{{*/
29 int Krigingx(double** ppredictions,double* obs_x, double* obs_y, double* obs_list, int obs_length,double* x_interp,double* y_interp,int n_interp,Options* options){
30
31 /*output*/
32 double *predictions = NULL;
33
34 /*Intermediaries*/
35- int i,j,n_obs;
36 double range;
37- double numerator,denominator,ratio;
38- double *x = NULL;
39- double *y = NULL;
40- double *obs = NULL;
41- double *Gamma = NULL;
42- double *GinvG0 = NULL;
43- double *Ginv1 = NULL;
44- double *GinvZ = NULL;
45- double *gamma0 = NULL;
46- double *ones = NULL;
47 char *output = NULL;
48 Variogram *variogram = NULL;
49 Observations *observations = NULL;
50
51+ /*threading: */
52+ KrigingxThreadStruct gate;
53+ int num=1;
54+#ifdef _MULTITHREADING_
55+ num=_NUMTHREADS_;
56+#endif
57+
58 /*Get Variogram from Options*/
59 ProcessVariogram(&variogram,options);
60 options->Get(&range,"searchrange",0.);
61@@ -44,29 +41,88 @@
62 /*Process observation dataset*/
63 observations=new Observations(obs_list,obs_x,obs_y,obs_length,options);
64
65- /*Allocation output*/
66+ /*Allocate output*/
67 predictions =(double*)xmalloc(n_interp*sizeof(double));
68- for(i=0;i<n_interp;i++) predictions[i]=0;
69+ for(int i=0;i<n_interp;i++) predictions[i]=0;
70
71 /*Get output*/
72- options->Get(&output,"output","quadtree");
73+ options->Get(&output,"output","prediction");
74
75 if(strcmp(output,"quadtree")==0){
76 observations->QuadtreeColoring(predictions,x_interp,y_interp,n_interp);
77- *ppredictions=predictions;
78- return 1;
79 }
80- if(strcmp(output,"variomap")==0){
81+ else if(strcmp(output,"variomap")==0){
82 observations->Variomap(predictions,x_interp,n_interp);
83- *ppredictions=predictions;
84- return 1;
85 }
86+ else if(strcmp(output,"prediction")==0){
87
88- /*Loop over all interpolations*/
89- printf(" interpolation progress: %5.2lf %%",0.0);
90- for(int idx=0;idx<n_interp;idx++){
91- if(idx%10==0) printf("\b\b\b\b\b\b\b%5.2lf %%",(double)idx/n_interp*100);
92+ /*initialize thread parameters: */
93+ gate.n_interp = n_interp;
94+ gate.x_interp = x_interp;
95+ gate.y_interp = y_interp;
96+ gate.range = range;
97+ gate.variogram = variogram;
98+ gate.observations = observations;
99+ gate.predictions = predictions;
100
101+ /*launch the thread manager with Krigingxt as a core: */
102+ LaunchThread(Krigingxt,(void*)&gate,num);
103+ }
104+ else{
105+ _error_("output '%s' not supported yet",output);
106+ }
107+
108+ /*clean-up and Assign output pointer*/
109+ delete variogram;
110+ delete observations;
111+ xfree((void**)&output);
112+ *ppredictions=predictions;
113+ return 1;
114+}/*}}}*/
115+/*FUNCTION Krigingxt{{{*/
116+void* Krigingxt(void* vpthread_handle){
117+
118+ /*gate variables :*/
119+ KrigingxThreadStruct *gate = NULL;
120+ pthread_handle *handle = NULL;
121+ int my_thread;
122+ int num_threads;
123+ int i0,i1;
124+
125+ /*recover handle and gate: */
126+ handle = (pthread_handle*)vpthread_handle;
127+ gate = (KrigingxThreadStruct*)handle->gate;
128+ my_thread = handle->id;
129+ num_threads = handle->num;
130+
131+ /*recover parameters :*/
132+ int n_interp = gate->n_interp;
133+ double *x_interp = gate->x_interp;
134+ double *y_interp = gate->y_interp;
135+ double range = gate->range;
136+ Variogram *variogram = gate->variogram;
137+ Observations *observations = gate->observations;
138+ double *predictions = gate->predictions;
139+
140+ /*Intermediaries*/
141+ int i,j,n_obs;
142+ double numerator,denominator,ratio;
143+ double *x = NULL;
144+ double *y = NULL;
145+ double *obs = NULL;
146+ double *Gamma = NULL;
147+ double *GinvG0 = NULL;
148+ double *Ginv1 = NULL;
149+ double *GinvZ = NULL;
150+ double *gamma0 = NULL;
151+ double *ones = NULL;
152+
153+ /*partition loop across threads: */
154+ if(my_thread==0) printf(" interpolation progress: %5.2lf %%",0.0);
155+ PartitionRange(&i0,&i1,n_interp,num_threads,my_thread);
156+ for(int idx=i0;idx<i1;idx++){
157+ if(my_thread==0 && idx%10==0) printf("\b\b\b\b\b\b\b%5.2lf %%",double(idx-i0)/double(i1-i0)*100);
158+
159 /*Get list of observations for current point*/
160 observations->ObservationList(&x,&y,&obs,&n_obs,x_interp[idx],y_interp[idx],range);
161
162@@ -88,9 +144,9 @@
163 for(i=0;i<n_obs;i++) gamma0[i] = variogram->SemiVariogram(x[i]-x_interp[idx],y[i]-y_interp[idx]);
164
165 /*Solve the three linear systems*/
166- GslSolve(&GinvG0,Gamma,gamma0,n_obs); // Gamma^-1 gamma0
167- GslSolve(&Ginv1, Gamma,ones,n_obs); // Gamma^-1 ones
168- GslSolve(&GinvZ, Gamma,obs,n_obs); // Gamma^-1 Z
169+ GslSolve(&GinvG0,Gamma,gamma0,n_obs); // Gamma^-1 gamma0
170+ GslSolve(&Ginv1, Gamma,ones,n_obs); // Gamma^-1 ones
171+ GslSolve(&GinvZ, Gamma,obs,n_obs); // Gamma^-1 Z
172
173 /*Prepare predictor*/
174 numerator=-1.; denominator=0.;
175@@ -112,15 +168,10 @@
176 xfree((void**)&Ginv1);
177 xfree((void**)&GinvZ);
178 }
179- printf("\b\b\b\b\b\b\b\b%5.2lf %%\n",100.0);
180+ if(my_thread==0) printf("\b\b\b\b\b\b\b\b%5.2lf %%\n",100.0);
181
182- /*clean-up and Assign output pointer*/
183- delete variogram;
184- delete observations;
185- xfree((void**)&output);
186- *ppredictions=predictions;
187- return 1;
188-}
189+ return NULL;
190+}/*}}}*/
191
192 void ProcessVariogram(Variogram **pvariogram,Options* options){/*{{{*/
193
194Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/modules/Krigingx/Krigingx.h
195===================================================================
196--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/modules/Krigingx/Krigingx.h (revision 12230)
197+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/modules/Krigingx/Krigingx.h (revision 12231)
198@@ -8,8 +8,23 @@
199 #include "../../objects/objects.h"
200 #include "../../toolkits/toolkits.h"
201
202+class Observations;
203+class Variogram;
204+
205 int Krigingx(double** ppredictions,double* x, double* y, double* observations, int n_obs,double* x_interp,double* y_interp,int n_interp,Options* options);
206 void ProcessVariogram(Variogram **pvariogram,Options* options);
207 void GslSolve(double** pX,double* A,double* B,int n);
208
209+/*threading: */
210+typedef struct{
211+ int n_interp;
212+ double *x_interp;
213+ double *y_interp;
214+ double range;
215+ Variogram *variogram;
216+ Observations *observations;
217+ double *predictions;
218+}KrigingxThreadStruct;
219+
220+void* Krigingxt(void*);
221 #endif /* _KRIGINGX_H */
222Index: /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/modules/Shp2Kmlx/Shp2Kmlx.h
223===================================================================
224--- /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/modules/Shp2Kmlx/Shp2Kmlx.h (revision 12230)
225+++ /proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/modules/Shp2Kmlx/Shp2Kmlx.h (revision 12231)
226@@ -21,9 +21,7 @@
227 #include "../../objects/objects.h"
228
229 /* local prototypes: */
230-int Shp2Kmlx(char* filshp,char* filkml,
231- int sgn);
232-int Shp2Kmlx(char* filshp,char* filkml,
233- int sgn,double cm,double sp);
234+int Shp2Kmlx(char* filshp,char* filkml, int sgn);
235+int Shp2Kmlx(char* filshp,char* filkml, int sgn,double cm,double sp);
236
237 #endif /* _SHP2KMLX_H */
Note: See TracBrowser for help on using the repository browser.