Ice Sheet System Model  4.18
Code documentation
Data Structures | Functions
Krigingx.h File Reference
#include "../../classes/classes.h"
#include "../../toolkits/toolkits.h"

Go to the source code of this file.

Data Structures

struct  KrigingxThreadStruct
 

Functions

int Krigingx (double **ppredictions, double **perror, double *x, double *y, double *observations, int n_obs, double *x_interp, double *y_interp, int n_interp, Options *options)
 
int pKrigingx (double **ppredictions, double **perror, double *x, double *y, double *observations, int n_obs, double *x_interp, double *y_interp, int n_interp, Options *options)
 
void ProcessVariogram (Variogram **pvariogram, Options *options)
 
void ProcessVariogram2 (Variogram **pvariogram, Options *options)
 
void * Krigingxt (void *)
 
void * NearestNeighbort (void *)
 
void * idwt (void *)
 
void * v4t (void *)
 
void * Distancest (void *)
 

Function Documentation

◆ Krigingx()

int Krigingx ( double **  ppredictions,
double **  perror,
double *  x,
double *  y,
double *  observations,
int  n_obs,
double *  x_interp,
double *  y_interp,
int  n_interp,
Options *  options 
)

Definition at line 13 of file Krigingx.cpp.

13  {/*{{{*/
14 
15  /*output*/
16  double *predictions = NULL;
17  double *error = NULL;
18 
19  /*Intermediaries*/
20  int mindata,maxdata;
21  double dmindata,dmaxdata,dnumthreads; //FIXME (Options come as double but we want to retrive integers)
22  double radius;
23  char *output = NULL;
24  Variogram *variogram = NULL;
25  Observations *observations = NULL;
26 
27  /*threading: */
28  KrigingxThreadStruct gate;
29  int num = _NUMTHREADS_;
30 
31  /*Get Variogram from Options*/
32  ProcessVariogram(&variogram,options);
33  options->Get(&radius,"searchradius",0.);
34  options->Get(&dmindata,"mindata",1.); mindata=int(dmindata);//FIXME (Options come as double but we want to retrive integers)
35  options->Get(&dmaxdata,"maxdata",50.); maxdata=int(dmaxdata);//FIXME (Options come as double but we want to retrive integers)
36  options->Get(&dnumthreads,"numthreads",double(num)); num=int(dnumthreads);//FIXME (Options come as double but we want to retrive integers)
37 
38  /*Process observation dataset*/
39  observations=new Observations(obs_list,obs_x,obs_y,obs_length,options);
40 
41  /*Allocate output*/
42  predictions =xNewZeroInit<double>(n_interp);
43  error =xNewZeroInit<double>(n_interp);
44 
45  /*Get output*/
46  options->Get(&output,"output",(char*)"prediction");
47 
48  if(strcmp(output,"quadtree")==0){
49  observations->QuadtreeColoring(predictions,x_interp,y_interp,n_interp);
50  }
51  else if(strcmp(output,"variomap")==0){
52  observations->Variomap(predictions,x_interp,n_interp);
53  }
54  else if(strcmp(output,"distance")==0){
55  /*initialize thread parameters: */
56  gate.n_interp = n_interp;
57  gate.x_interp = x_interp;
58  gate.y_interp = y_interp;
59  gate.radius = radius;
60  gate.mindata = mindata;
61  gate.maxdata = maxdata;
62  gate.variogram = variogram;
63  gate.observations = observations;
64  gate.predictions = predictions;
65  gate.error = error;
66  gate.numdone = xNewZeroInit<int>(num);
67 
68  /*launch the thread manager with Krigingxt as a core: */
69  LaunchThread(Distancest,(void*)&gate,num);
70  xDelete<int>(gate.numdone);
71  }
72  else if(strcmp(output,"delaunay")==0){
73 
74  #ifdef _HAVE_BAMG_
75  int nobs,nel;
76  double *x = NULL;
77  double *y = NULL;
78  double *data = NULL;
79  int *index = NULL;
80 
81  observations->ObservationList(&x,&y,&data,&nobs);
82 
83  _printf_("Generation Delaunay Triangulation\n");
84  BamgTriangulatex(&index,&nel,x,y,nobs);
85 
86  _printf_("Interpolating\n");
87  xDelete<double>(predictions);
88  InterpFromMeshToMesh2dx(&predictions,index,x,y,nobs,nel,data,nobs,1,x_interp,y_interp,n_interp,options);
89  xDelete<double>(x);
90  xDelete<double>(y);
91  xDelete<double>(data);
92  xDelete<int>(index);
93  #else
94  _error_("you did not compile ISSM with bamg");
95  #endif
96  }
97  else if(strcmp(output,"nearestneighbor")==0){
98  /*initialize thread parameters: */
99  gate.n_interp = n_interp;
100  gate.x_interp = x_interp;
101  gate.y_interp = y_interp;
102  gate.radius = radius;
103  gate.mindata = mindata;
104  gate.maxdata = maxdata;
105  gate.variogram = variogram;
106  gate.observations = observations;
107  gate.predictions = predictions;
108  gate.error = error;
109  gate.numdone = xNewZeroInit<int>(num);
110 
111  /*launch the thread manager with Krigingxt as a core: */
112  LaunchThread(NearestNeighbort,(void*)&gate,num);
113  _printf_("\r interpolation progress: "<<fixed<<setw(6)<<setprecision(2)<<100.<<"% \n");
114  xDelete<int>(gate.numdone);
115  }
116  else if(strcmp(output,"idw")==0){ //Inverse distance weighting
117  double power;
118  options->Get(&power,"power",2.);
119  /*initialize thread parameters: */
120  gate.n_interp = n_interp;
121  gate.x_interp = x_interp;
122  gate.y_interp = y_interp;
123  gate.radius = radius;
124  gate.mindata = mindata;
125  gate.maxdata = maxdata;
126  gate.variogram = variogram;
127  gate.observations = observations;
128  gate.predictions = predictions;
129  gate.error = error;
130  gate.numdone = xNewZeroInit<int>(num);
131  gate.power = power;
132 
133  /*launch the thread manager with Krigingxt as a core: */
134  LaunchThread(idwt,(void*)&gate,num);
135  _printf_("\r interpolation progress: "<<fixed<<setw(6)<<setprecision(2)<<100.<<"% \n");
136  xDelete<int>(gate.numdone);
137  }
138  else if(strcmp(output,"v4")==0){ //Inverse distance weighting
139 #if !defined(_HAVE_GSL_)
140  _error_("GSL is required for v4 interpolation");
141 #endif
142  /*initialize thread parameters: */
143  gate.n_interp = n_interp;
144  gate.x_interp = x_interp;
145  gate.y_interp = y_interp;
146  gate.radius = radius;
147  gate.mindata = mindata;
148  gate.maxdata = maxdata;
149  gate.variogram = variogram;
150  gate.observations = observations;
151  gate.predictions = predictions;
152  gate.error = error;
153  gate.numdone = xNewZeroInit<int>(num);
154 
155  /*launch the thread manager with Krigingxt as a core: */
156  LaunchThread(v4t,(void*)&gate,num);
157  _printf_("\r interpolation progress: "<<fixed<<setw(6)<<setprecision(2)<<100.<<"% \n");
158  xDelete<int>(gate.numdone);
159  }
160  else if(strcmp(output,"prediction")==0){
161 #if !defined(_HAVE_GSL_)
162  _error_("GSL is required for v4 interpolation");
163 #endif
164 
165  /*initialize thread parameters: */
166  gate.n_interp = n_interp;
167  gate.x_interp = x_interp;
168  gate.y_interp = y_interp;
169  gate.radius = radius;
170  gate.mindata = mindata;
171  gate.maxdata = maxdata;
172  gate.variogram = variogram;
173  gate.observations = observations;
174  gate.predictions = predictions;
175  gate.error = error;
176  gate.numdone = xNewZeroInit<int>(num);
177 
178  /*launch the thread manager with Krigingxt as a core: */
179  LaunchThread(Krigingxt,(void*)&gate,num);
180  _printf_("\r interpolation progress: "<<fixed<<setw(6)<<setprecision(2)<<100.<<"% \n");
181  xDelete<int>(gate.numdone);
182  }
183  else{
184  _error_("output '" << output << "' not supported yet");
185  }
186 
187  /*clean-up and Assign output pointer*/
188  delete variogram;
189  delete observations;
190  xDelete<char>(output);
191  *ppredictions = predictions;
192  *perror = error;
193  return 1;
194 }/*}}}*/

◆ pKrigingx()

int pKrigingx ( double **  ppredictions,
double **  perror,
double *  x,
double *  y,
double *  observations,
int  n_obs,
double *  x_interp,
double *  y_interp,
int  n_interp,
Options *  options 
)

Definition at line 11 of file pKrigingx.cpp.

11  {/*{{{*/
12 
13 #ifdef _HAVE_MPI_
14  int num_procs;
15  int my_rank;
16 
17  /*output*/
18  double *predictions = NULL;
19  double *error = NULL;
20 
21  /*Intermediaries*/
22  int mindata,maxdata;
23  double dmindata,dmaxdata;
24  double radius;
25  char *output = NULL;
26  Variogram *variogram = NULL;
27  Observations *observations = NULL;
28 
29  /*timing*/
30  double start, finish;
31  double start_core, finish_core;
32  double start_init, finish_init;
33 
34  /*Get my_rank: */
35  my_rank=IssmComm::GetRank();
36  num_procs=IssmComm::GetSize();
37 
38  /*Get some Options*/
39  ISSM_MPI_Barrier(ISSM_MPI_COMM_WORLD); start=ISSM_MPI_Wtime();
40  options->Get(&radius,"searchradius",0.);
41 
42  options->Get(&dmindata,"mindata",1.); mindata=int(dmindata);//FIXME (Options come as double but we want to retrive integers)
43  options->Get(&dmaxdata,"maxdata",50.); maxdata=int(dmaxdata);//FIXME (Options come as double but we want to retrive integers)
44 
45  /*Process observation dataset*/
46  ISSM_MPI_Barrier(ISSM_MPI_COMM_WORLD); start_init=ISSM_MPI_Wtime();
47  observations=new Observations(obs_list,obs_x,obs_y,obs_length,options);
48  ISSM_MPI_Barrier(ISSM_MPI_COMM_WORLD); finish_init=ISSM_MPI_Wtime();
49 
50  /*Allocate output*/
51  predictions =xNewZeroInit<double>(n_interp);
52  error =xNewZeroInit<double>(n_interp);
53 
54  /*Get output*/
55  options->Get(&output,"output",(char*)"prediction");
56 
57  ISSM_MPI_Barrier(ISSM_MPI_COMM_WORLD); start_core=ISSM_MPI_Wtime( );
58  if(strcmp(output,"quadtree")==0){
59  observations->QuadtreeColoring(predictions,x_interp,y_interp,n_interp);
60  }
61  else if(strcmp(output,"variomap")==0){
62  observations->Variomap(predictions,x_interp,n_interp);
63  }
64  else if(strcmp(output,"prediction")==0){
65 
66  /*Process Variogram*/
67  ProcessVariogram2(&variogram,options);
68 
69  /*partition loop across threads: */
70  for(int idx=my_rank;idx<n_interp;idx+=num_procs){
71  _printf0_(" interpolation progress: "<<fixed<<setw(6)<<setprecision(4)<<double(idx)/double(n_interp)*100.<<"% \n");
72  observations->InterpolationKriging(&predictions[idx],&error[idx],x_interp[idx],y_interp[idx],radius,mindata,maxdata,variogram);
73  }
74  _printf0_(" interpolation progress: "<<fixed<<setw(6)<<setprecision(4)<<100.<<"% \n");
75 
76  double *sumpredictions =xNew<double>(n_interp);
77  double *sumerror =xNew<double>(n_interp);
78  ISSM_MPI_Allreduce(predictions,sumpredictions,n_interp,ISSM_MPI_PDOUBLE,ISSM_MPI_SUM,IssmComm::GetComm());
79  ISSM_MPI_Allreduce(error,sumerror,n_interp,ISSM_MPI_PDOUBLE,ISSM_MPI_SUM,IssmComm::GetComm());
80  xDelete<double>(error); error=sumerror;
81  xDelete<double>(predictions); predictions=sumpredictions;
82  }
83  else if(strcmp(output,"v4")==0){
84 
85  /*partition loop across threads: */
86  for(int idx=my_rank;idx<n_interp;idx+=num_procs){
87  _printf0_(" interpolation progress: "<<fixed<<setw(6)<<setprecision(4)<<double(idx)/double(n_interp)*100.<<"% \n");
88  observations->InterpolationV4(&predictions[idx],x_interp[idx],y_interp[idx],radius,mindata,maxdata);
89  }
90  _printf0_(" interpolation progress: "<<fixed<<setw(6)<<setprecision(4)<<100.<<"% \n");
91 
92  double *sumpredictions =xNew<double>(n_interp);
93  ISSM_MPI_Allreduce(predictions,sumpredictions,n_interp,ISSM_MPI_PDOUBLE,ISSM_MPI_SUM,IssmComm::GetComm());
94  xDelete<double>(predictions); predictions=sumpredictions;
95  }
96  else if(strcmp(output,"nearestneighbor")==0){
97 
98  /*partition loop across threads: */
99  for(int idx=my_rank;idx<n_interp;idx+=num_procs){
100  _printf0_(" interpolation progress: "<<fixed<<setw(6)<<setprecision(4)<<double(idx)/double(n_interp)*100.<<"% \n");
101  observations->InterpolationNearestNeighbor(&predictions[idx],x_interp[idx],y_interp[idx],radius);
102  }
103  _printf0_(" interpolation progress: "<<fixed<<setw(6)<<setprecision(4)<<100.<<"% \n");
104 
105  double *sumpredictions =xNew<double>(n_interp);
106  ISSM_MPI_Allreduce(predictions,sumpredictions,n_interp,ISSM_MPI_PDOUBLE,ISSM_MPI_SUM,IssmComm::GetComm());
107  xDelete<double>(predictions); predictions=sumpredictions;
108  }
109  else if(strcmp(output,"distance")==0){
110 
111  /*partition loop across threads: */
112  for(int idx=my_rank;idx<n_interp;idx+=num_procs){
113  _printf0_(" interpolation progress: "<<fixed<<setw(6)<<setprecision(4)<<double(idx)/double(n_interp)*100.<<"% \n");
114  observations->Distances(&predictions[idx],&x_interp[idx],&y_interp[idx],1,radius);
115  }
116  _printf0_(" interpolation progress: "<<fixed<<setw(6)<<setprecision(4)<<100.<<"% \n");
117 
118  double *sumpredictions =xNew<double>(n_interp);
119  ISSM_MPI_Allreduce(predictions,sumpredictions,n_interp,ISSM_MPI_PDOUBLE,ISSM_MPI_SUM,IssmComm::GetComm());
120  xDelete<double>(predictions); predictions=sumpredictions;
121  }
122  else if(strcmp(output,"idw")==0){
123  double power;
124  options->Get(&power,"power",2.);
125 
126  /*partition loop across threads: */
127  for(int idx=my_rank;idx<n_interp;idx+=num_procs){
128  _printf0_(" interpolation progress: "<<fixed<<setw(6)<<setprecision(4)<<double(idx)/double(n_interp)*100.<<"% \n");
129  observations->InterpolationIDW(&predictions[idx],x_interp[idx],y_interp[idx],radius,mindata,maxdata,power);
130  }
131  _printf0_(" interpolation progress: "<<fixed<<setw(6)<<setprecision(4)<<100.<<"% \n");
132 
133  double *sumpredictions =xNew<double>(n_interp);
134  ISSM_MPI_Allreduce(predictions,sumpredictions,n_interp,ISSM_MPI_PDOUBLE,ISSM_MPI_SUM,IssmComm::GetComm());
135  xDelete<double>(predictions); predictions=sumpredictions;
136  }
137  else{
138  _error_("output '" << output << "' not supported yet");
139  }
140  ISSM_MPI_Barrier(ISSM_MPI_COMM_WORLD); finish_core=ISSM_MPI_Wtime( );
141 
142  /*clean-up and Assign output pointer*/
143  delete variogram;
144  delete observations;
145  xDelete<char>(output);
146  *ppredictions = predictions;
147  *perror = error;
148 
149  ISSM_MPI_Barrier(ISSM_MPI_COMM_WORLD); finish=ISSM_MPI_Wtime( );
150  _printf0_("\n " << setw(34) << left << "Observation fitering elapsed time: " << finish_init-start_init << " seconds \n\n");
151  _printf0_(" " << setw(34) << left << "Kriging prediction elapsed time: " << finish_core-start_core << " seconds \n\n");
152  _printf0_("\n " << "Total elapsed time " << int((finish-start)/3600) << " hrs " << int(int(finish-start)%3600/60) << " min " << int(finish-start)%60 << " sec\n\n\n");
153  return 1;
154 #else
155  _error_("MPI not available");
156 #endif
157 }/*}}}*/

◆ ProcessVariogram()

void ProcessVariogram ( Variogram **  pvariogram,
Options *  options 
)

Definition at line 409 of file Krigingx.cpp.

409  {/*{{{*/
410 
411  /*Intermediaries*/
412  Variogram* variogram = NULL;
413  char *model = NULL;
414 
415  if(options->GetOption("model")){
416  options->Get(&model,"model");
417  if (strcmp(model,"gaussian")==0) variogram = new GaussianVariogram(options);
418  else if(strcmp(model,"exponential")==0) variogram = new ExponentialVariogram(options);
419  else if(strcmp(model,"spherical")==0) variogram = new SphericalVariogram(options);
420  else if(strcmp(model,"power")==0) variogram = new PowerVariogram(options);
421  else _error_("variogram " << model << " not supported yet (list of supported variogram: gaussian, exponential, spherical and power)");
422  }
423  else variogram = new GaussianVariogram(options);
424 
425  /*Assign output pointer*/
426  xDelete<char>(model);
427  *pvariogram = variogram;
428 }/*}}}*/

◆ ProcessVariogram2()

void ProcessVariogram2 ( Variogram **  pvariogram,
Options *  options 
)

Definition at line 158 of file pKrigingx.cpp.

158  {/*{{{*/
159 
160  /*Intermediaries*/
161  Variogram* variogram = NULL;
162  char *model = NULL;
163 
164  if(options->GetOption("model")){
165  options->Get(&model,"model");
166  if (strcmp(model,"gaussian")==0) variogram = new GaussianVariogram(options);
167  else if(strcmp(model,"exponential")==0) variogram = new ExponentialVariogram(options);
168  else if(strcmp(model,"spherical")==0) variogram = new SphericalVariogram(options);
169  else if(strcmp(model,"power")==0) variogram = new PowerVariogram(options);
170  else _error_("variogram " << model << " not supported yet (list of supported variogram: gaussian, exponential, spherical and power)");
171  }
172  else variogram = new GaussianVariogram(options);
173 
174  /*Assign output pointer*/
175  xDelete<char>(model);
176  *pvariogram = variogram;
177 }/*}}}*/

◆ Krigingxt()

void* Krigingxt ( void *  )

Definition at line 195 of file Krigingx.cpp.

195  {/*{{{*/
196 
197  /*gate variables :*/
198  KrigingxThreadStruct *gate = NULL;
199  pthread_handle *handle = NULL;
200  int my_thread;
201  int num_threads;
202  int i0,i1;
203 
204  /*recover handle and gate: */
205  handle = (pthread_handle*)vpthread_handle;
206  gate = (KrigingxThreadStruct*)handle->gate;
207  my_thread = handle->id;
208  num_threads = handle->num;
209 
210  /*recover parameters :*/
211  int n_interp = gate->n_interp;
212  double *x_interp = gate->x_interp;
213  double *y_interp = gate->y_interp;
214  double radius = gate->radius;
215  int mindata = gate->mindata;
216  int maxdata = gate->maxdata;
217  Variogram *variogram = gate->variogram;
218  Observations *observations = gate->observations;
219  double *predictions = gate->predictions;
220  double *error = gate->error;
221  int *numdone = gate->numdone;
222 
223  /*partition loop across threads: */
224  PartitionRange(&i0,&i1,n_interp,num_threads,my_thread);
225  for(int idx=i0;idx<i1;idx++){
226 
227  /*Print info*/
228  numdone[my_thread]=idx-i0;
229  if(my_thread==0){
230  int alldone=numdone[0];
231  for(int i=1;i<num_threads;i++) alldone+=numdone[i];
232  _printf_("\r interpolation progress: "<<setw(6)<<setprecision(2)<<double(alldone)/double(n_interp)*100.<<"% ");
233  }
234 
235  /*Kriging interpolation*/
236  observations->InterpolationKriging(&predictions[idx],&error[idx],x_interp[idx],y_interp[idx],radius,mindata,maxdata,variogram);
237  }
238 
239  return NULL;
240 }/*}}}*/

◆ NearestNeighbort()

void* NearestNeighbort ( void *  )

Definition at line 241 of file Krigingx.cpp.

241  {/*{{{*/
242 
243  /*gate variables :*/
244  KrigingxThreadStruct *gate = NULL;
245  pthread_handle *handle = NULL;
246  int my_thread;
247  int num_threads;
248  int i0,i1;
249 
250  /*recover handle and gate: */
251  handle = (pthread_handle*)vpthread_handle;
252  gate = (KrigingxThreadStruct*)handle->gate;
253  my_thread = handle->id;
254  num_threads = handle->num;
255 
256  /*recover parameters :*/
257  int n_interp = gate->n_interp;
258  double *x_interp = gate->x_interp;
259  double *y_interp = gate->y_interp;
260  double radius = gate->radius;
261  int mindata = gate->mindata;
262  int maxdata = gate->maxdata;
263  Variogram *variogram = gate->variogram;
264  Observations *observations = gate->observations;
265  double *predictions = gate->predictions;
266  double *error = gate->error;
267  int *numdone = gate->numdone;
268 
269  /*partition loop across threads: */
270  PartitionRange(&i0,&i1,n_interp,num_threads,my_thread);
271  for(int idx=i0;idx<i1;idx++){
272 
273  /*Print info*/
274  numdone[my_thread]=idx-i0;
275  if(my_thread==0){
276  int alldone=numdone[0];
277  for(int i=1;i<num_threads;i++) alldone+=numdone[i];
278  _printf_("\r interpolation progress: "<<setw(6)<<setprecision(2)<<double(alldone)/double(n_interp)*100.<<"% ");
279  }
280 
281  observations->InterpolationNearestNeighbor(&predictions[idx],x_interp[idx],y_interp[idx],radius);
282  }
283 
284  return NULL;
285 }/*}}}*/

◆ idwt()

void* idwt ( void *  )

Definition at line 286 of file Krigingx.cpp.

286  {/*{{{*/
287 
288  /*gate variables :*/
289  KrigingxThreadStruct *gate = NULL;
290  pthread_handle *handle = NULL;
291  int my_thread;
292  int num_threads;
293  int i0,i1;
294 
295  /*recover handle and gate: */
296  handle = (pthread_handle*)vpthread_handle;
297  gate = (KrigingxThreadStruct*)handle->gate;
298  my_thread = handle->id;
299  num_threads = handle->num;
300 
301  /*recover parameters :*/
302  int n_interp = gate->n_interp;
303  double *x_interp = gate->x_interp;
304  double *y_interp = gate->y_interp;
305  double radius = gate->radius;
306  int mindata = gate->mindata;
307  int maxdata = gate->maxdata;
308  Variogram *variogram = gate->variogram;
309  Observations *observations = gate->observations;
310  double *predictions = gate->predictions;
311  double *error = gate->error;
312  int *numdone = gate->numdone;
313  double power = gate->power;
314 
315  /*partition loop across threads: */
316  PartitionRange(&i0,&i1,n_interp,num_threads,my_thread);
317  for(int idx=i0;idx<i1;idx++){
318 
319  /*Print info*/
320  numdone[my_thread]=idx-i0;
321  if(my_thread==0){
322  int alldone=numdone[0];
323  for(int i=1;i<num_threads;i++) alldone+=numdone[i];
324  _printf_("\r interpolation progress: "<<setw(6)<<setprecision(2)<<double(alldone)/double(n_interp)*100.<<"% ");
325  }
326 
327  observations->InterpolationIDW(&predictions[idx],x_interp[idx],y_interp[idx],radius,mindata,maxdata,power);
328  }
329  return NULL;
330 }/*}}}*/

◆ v4t()

void* v4t ( void *  )

Definition at line 331 of file Krigingx.cpp.

331  {/*{{{*/
332 
333  /*gate variables :*/
334  KrigingxThreadStruct *gate = NULL;
335  pthread_handle *handle = NULL;
336  int my_thread;
337  int num_threads;
338  int i0,i1;
339 
340  /*recover handle and gate: */
341  handle = (pthread_handle*)vpthread_handle;
342  gate = (KrigingxThreadStruct*)handle->gate;
343  my_thread = handle->id;
344  num_threads = handle->num;
345 
346  /*recover parameters :*/
347  int n_interp = gate->n_interp;
348  double *x_interp = gate->x_interp;
349  double *y_interp = gate->y_interp;
350  double radius = gate->radius;
351  int mindata = gate->mindata;
352  int maxdata = gate->maxdata;
353  Variogram *variogram = gate->variogram;
354  Observations *observations = gate->observations;
355  double *predictions = gate->predictions;
356  double *error = gate->error;
357  int *numdone = gate->numdone;
358 
359  /*partition loop across threads: */
360  PartitionRange(&i0,&i1,n_interp,num_threads,my_thread);
361  for(int idx=i0;idx<i1;idx++){
362 
363  /*Print info*/
364  numdone[my_thread]=idx-i0;
365  if(my_thread==0){
366  int alldone=numdone[0];
367  for(int i=1;i<num_threads;i++) alldone+=numdone[i];
368  _printf_("\r interpolation progress: "<<setw(6)<<setprecision(2)<<double(alldone)/double(n_interp)*100.<<"% ");
369  }
370 
371  observations->InterpolationV4(&predictions[idx],x_interp[idx],y_interp[idx],radius,mindata,maxdata);
372  }
373  return NULL;
374 }/*}}}*/

◆ Distancest()

void* Distancest ( void *  )

Definition at line 375 of file Krigingx.cpp.

375  {/*{{{*/
376 
377  /*gate variables :*/
378  KrigingxThreadStruct *gate = NULL;
379  pthread_handle *handle = NULL;
380  int my_thread;
381  int num_threads;
382  int i0,i1;
383 
384  /*recover handle and gate: */
385  handle = (pthread_handle*)vpthread_handle;
386  gate = (KrigingxThreadStruct*)handle->gate;
387  my_thread = handle->id;
388  num_threads = handle->num;
389 
390  /*recover parameters :*/
391  int n_interp = gate->n_interp;
392  double *x_interp = gate->x_interp;
393  double *y_interp = gate->y_interp;
394  double radius = gate->radius;
395  int mindata = gate->mindata;
396  int maxdata = gate->maxdata;
397  Variogram *variogram = gate->variogram;
398  Observations *observations = gate->observations;
399  double *predictions = gate->predictions;
400  double *error = gate->error;
401  int *numdone = gate->numdone;
402 
403  /*partition loop across threads: */
404  PartitionRange(&i0,&i1,n_interp,num_threads,my_thread);
405  observations->Distances(&predictions[i0],&x_interp[i0],&y_interp[i0],i1-i0,radius);
406  return NULL;
407 }/*}}}*/
KrigingxThreadStruct::observations
Observations * observations
Definition: Krigingx.h:28
KrigingxThreadStruct::variogram
Variogram * variogram
Definition: Krigingx.h:27
Observations::Distances
void Distances(IssmPDouble *distances, IssmPDouble *x, IssmPDouble *y, int n, IssmPDouble radius)
Definition: Observations.cpp:302
pthread_handle::id
int id
Definition: issm_threads.h:12
_printf0_
#define _printf0_(StreamArgs)
Definition: Print.h:29
ISSM_MPI_Allreduce
int ISSM_MPI_Allreduce(void *sendbuf, void *recvbuf, int count, ISSM_MPI_Datatype datatype, ISSM_MPI_Op op, ISSM_MPI_Comm comm)
Definition: issmmpi.cpp:116
v4t
void * v4t(void *vpthread_handle)
Definition: Krigingx.cpp:331
_printf_
#define _printf_(StreamArgs)
Definition: Print.h:22
LaunchThread
void LaunchThread(void *function(void *), void *gate, int num_threads)
Definition: LaunchThread.cpp:25
ISSM_MPI_SUM
#define ISSM_MPI_SUM
Definition: issmmpi.h:134
IssmComm::GetComm
static ISSM_MPI_Comm GetComm(void)
Definition: IssmComm.cpp:30
pthread_handle::num
int num
Definition: issm_threads.h:13
GaussianVariogram
Definition: GaussianVariogram.h:12
KrigingxThreadStruct
Definition: Krigingx.h:20
KrigingxThreadStruct::error
double * error
Definition: Krigingx.h:30
KrigingxThreadStruct::mindata
int mindata
Definition: Krigingx.h:25
KrigingxThreadStruct::y_interp
double * y_interp
Definition: Krigingx.h:23
KrigingxThreadStruct::n_interp
int n_interp
Definition: Krigingx.h:21
ProcessVariogram2
void ProcessVariogram2(Variogram **pvariogram, Options *options)
Definition: pKrigingx.cpp:158
ISSM_MPI_COMM_WORLD
#define ISSM_MPI_COMM_WORLD
Definition: issmmpi.h:137
Observations::InterpolationIDW
void InterpolationIDW(IssmDouble *pprediction, IssmDouble x_interp, IssmDouble y_interp, IssmDouble radius, int mindata, int maxdata, IssmDouble power)
Definition: Observations.cpp:481
PowerVariogram
Definition: PowerVariogram.h:11
Observations::QuadtreeColoring
void QuadtreeColoring(IssmDouble *A, IssmDouble *x, IssmDouble *y, int n)
Definition: Observations.cpp:692
BamgTriangulatex
int BamgTriangulatex(int **pindex, int *pnels, double *x, double *y, int nods)
Definition: BamgTriangulatex.cpp:14
Observations::Variomap
void Variomap(IssmDouble *gamma, IssmDouble *x, int n)
Definition: Observations.cpp:704
pthread_handle
Definition: issm_threads.h:10
IssmComm::GetSize
static int GetSize(void)
Definition: IssmComm.cpp:46
PartitionRange
void PartitionRange(int *pi0, int *pi1, int num_el, int num_threads, int my_thread)
Definition: PartitionRange.cpp:13
KrigingxThreadStruct::numdone
int * numdone
Definition: Krigingx.h:31
Krigingxt
void * Krigingxt(void *vpthread_handle)
Definition: Krigingx.cpp:195
ISSM_MPI_PDOUBLE
#define ISSM_MPI_PDOUBLE
Definition: issmmpi.h:126
Options::Get
void Get(OptionType *pvalue, const char *name)
Definition: Options.h:21
KrigingxThreadStruct::radius
double radius
Definition: Krigingx.h:24
Observations::ObservationList
void ObservationList(IssmDouble **px, IssmDouble **py, IssmDouble **pobs, int *pnobs)
Definition: Observations.cpp:316
IssmComm::GetRank
static int GetRank(void)
Definition: IssmComm.cpp:34
ProcessVariogram
void ProcessVariogram(Variogram **pvariogram, Options *options)
Definition: Krigingx.cpp:409
idwt
void * idwt(void *vpthread_handle)
Definition: Krigingx.cpp:286
Distancest
void * Distancest(void *vpthread_handle)
Definition: Krigingx.cpp:375
Options::GetOption
Option * GetOption(const char *name)
Definition: Options.cpp:67
_error_
#define _error_(StreamArgs)
Definition: exceptions.h:49
KrigingxThreadStruct::predictions
double * predictions
Definition: Krigingx.h:29
Observations::InterpolationKriging
void InterpolationKriging(IssmDouble *pprediction, IssmDouble *perror, IssmDouble x_interp, IssmDouble y_interp, IssmDouble radius, int mindata, int maxdata, Variogram *variogram)
Definition: Observations.cpp:526
NearestNeighbort
void * NearestNeighbort(void *vpthread_handle)
Definition: Krigingx.cpp:241
SphericalVariogram
Definition: SphericalVariogram.h:11
ISSM_MPI_Wtime
double ISSM_MPI_Wtime(void)
Definition: issmmpi.cpp:511
Variogram
Definition: Variogram.h:10
ExponentialVariogram
Definition: ExponentialVariogram.h:11
ISSM_MPI_Barrier
int ISSM_MPI_Barrier(ISSM_MPI_Comm comm)
Definition: issmmpi.cpp:148
KrigingxThreadStruct::power
double power
Definition: Krigingx.h:32
pthread_handle::gate
void * gate
Definition: issm_threads.h:11
Observations::InterpolationNearestNeighbor
void InterpolationNearestNeighbor(IssmDouble *pprediction, IssmDouble x_interp, IssmDouble y_interp, IssmDouble radius)
Definition: Observations.cpp:599
Observations::InterpolationV4
void InterpolationV4(IssmDouble *pprediction, IssmDouble x_interp, IssmDouble y_interp, IssmDouble radius, int mindata, int maxdata)
Definition: Observations.cpp:610
Observations
Declaration of Observations class.
Definition: Observations.h:16
KrigingxThreadStruct::x_interp
double * x_interp
Definition: Krigingx.h:22
InterpFromMeshToMesh2dx
int InterpFromMeshToMesh2dx(double **pdata_interp, int *index_data, double *x_data, double *y_data, int nods_data, int nels_data, double *data, int M_data, int N_data, double *x_interp, double *y_interp, int N_interp, Options *options)
Definition: InterpFromMeshToMesh2dx.cpp:16
KrigingxThreadStruct::maxdata
int maxdata
Definition: Krigingx.h:26
Generated on Thu Jul 2 2020 08:09:24 for Ice Sheet System Model by   doxygen 1.8.19