YetAnotherCoupler 3.5.2
Loading...
Searching...
No Matches
read_mpiom_grid.c
Go to the documentation of this file.
1// Copyright (c) 2024 The YAC Authors
2//
3// SPDX-License-Identifier: BSD-3-Clause
4
5#ifdef HAVE_CONFIG_H
6#include "config.h"
7#endif
8
9#include <stdlib.h>
10#include <stdio.h>
11#include <stdint.h>
12
13#include "read_mpiom_grid.h"
14#include "geometry.h"
15#include "utils_common.h"
16#include "io_utils.h"
17
18#ifdef YAC_NETCDF_ENABLED
19
20#include <netcdf.h>
21
22
23static void get_mpiom_vertices(int ncid,
24 double **vertex_lon, // longitude of vertex coordinates
25 double **vertex_lat, // latitude of vertex coordinates
26 size_t *nbr_vertices);
27
28static void get_mpiom_cell_center(int ncid,
29 double **cell_lon, // longitude of cell coordinates
30 double **cell_lat, // latitude of cell coordinates
31 size_t *nbr_cells);
32
33static void get_mpiom_cell_mask(int ncid,
34 int **cell_mask, // data for mask
35 size_t *nbr_cells );
36
37static void remove_duplicated_coords(double * temp_coords_lon,
38 double * temp_coords_lat,
39 int * temp_coords_mask,
40 size_t * temp_nbr_coords,
41 int * old_to_new_id);
42
43void yac_read_mpiom_grid_information(const char * filename, int * nbr_vertices,
44 int * nbr_cells, int ** num_vertices_per_cell,
45 int ** cell_to_vertex, double ** x_vertices,
46 double ** y_vertices, double ** x_cells,
47 double ** y_cells, int ** cell_mask) {
48
49 int ncid;
50
51 /* Open file */
52
53 yac_nc_open(filename, NC_NOWRITE, &ncid);
54
55 /* Get vertex longitudes and latitudes of cells and
56 relations between vertices and cells*/
57
58 {
59 double * temp_vertex_lon;
60 double * temp_vertex_lat;
61 size_t temp_nbr_vertices;
62
64 ncid, &temp_vertex_lon, &temp_vertex_lat, &temp_nbr_vertices);
65
66 *nbr_cells = (int)(temp_nbr_vertices / 4);
67
68 int * old_to_new_id = xmalloc(temp_nbr_vertices * sizeof(*old_to_new_id));
69
71 temp_vertex_lon, temp_vertex_lat, NULL,
72 &temp_nbr_vertices, old_to_new_id);
73
74 *nbr_vertices = temp_nbr_vertices;
75
76 *x_vertices =
77 xrealloc(temp_vertex_lon, temp_nbr_vertices * sizeof(*temp_vertex_lon));
78 *y_vertices =
79 xrealloc(temp_vertex_lat, temp_nbr_vertices * sizeof(*temp_vertex_lat));
80
81 for (size_t i = 0; i < temp_nbr_vertices; ++i) {
82 (*x_vertices)[i] *= YAC_RAD;
83 (*y_vertices)[i] *= YAC_RAD;
84 }
85
86 *num_vertices_per_cell = xmalloc(*nbr_cells * sizeof(**num_vertices_per_cell));
87 for (int i = 0; i < *nbr_cells; ++i) (*num_vertices_per_cell)[i] = 4;
88
89 *cell_to_vertex = xmalloc(*nbr_cells * 4 * sizeof(**cell_to_vertex));
90
91 // Unfortunately the data is only available in Fortran order
92 for (int i = 0; i < *nbr_cells; ++i)
93 for (int j = 0; j < 4; ++j)
94 (*cell_to_vertex)[4*i+j] = old_to_new_id[4*i+j];
95 free(old_to_new_id);
96 }
97
98 {
99 /* Get longitudes and latitudes of cell center points */
100
101 double * temp_cell_lon;
102 double * temp_cell_lat;
103 int * temp_cell_mask;
104
105 size_t temp_nbr_cells;
106 size_t temp_nbr_masks;
107
109 ncid, &temp_cell_lon, &temp_cell_lat, &temp_nbr_cells);
110
111 get_mpiom_cell_mask(ncid, &temp_cell_mask, &temp_nbr_masks);
112
114 temp_nbr_cells == temp_nbr_masks,
115 "ERROR(yac_read_mpiom_grid_information): "
116 "missmatch between cell count and mask value count")
117
118 int * old_to_new_id = xmalloc(temp_nbr_cells * sizeof(*old_to_new_id));
119
121 temp_cell_lon, temp_cell_lat, temp_cell_mask,
122 &temp_nbr_cells, old_to_new_id);
123
124 *nbr_cells = temp_nbr_cells;
125
126 *x_cells = xrealloc(temp_cell_lon, temp_nbr_cells * sizeof(**x_cells));
127 *y_cells = xrealloc(temp_cell_lat, temp_nbr_cells * sizeof(**y_cells));
128
129 for (size_t i = 0; i < temp_nbr_cells; ++i) {
130 (*x_cells)[i] *= YAC_RAD;
131 (*y_cells)[i] *= YAC_RAD;
132 }
133
134 *cell_mask = xrealloc(temp_cell_mask,temp_nbr_cells * sizeof(int));
135 for (size_t i = 0; i < temp_nbr_cells; ++i)
136 (*cell_mask)[i] = !(*cell_mask)[i];
137
138 free(old_to_new_id);
139 }
140 YAC_HANDLE_ERROR(nc_close(ncid));
141}
142
144 const char * filename, int * num_vertices, int * num_cells,
145 int ** num_vertices_per_cell, int ** cell_to_vertex, double ** x_vertices,
146 double ** y_vertices, double ** x_cells,
147 double ** y_cells, int ** global_cell_id,
148 int ** cell_mask, int ** cell_core_mask,
149 int ** global_corner_id, int ** corner_core_mask, int rank, int size) {
150
151 // read the global grid
152 yac_read_mpiom_grid_information(filename, num_vertices, num_cells,
153 num_vertices_per_cell, cell_to_vertex,
154 x_vertices, y_vertices, x_cells, y_cells, cell_mask);
155
156 // determine local rank
157 int num_cells_per_process = (*num_cells + size - 1)/size;
158 int local_start = rank * num_cells_per_process;
159 int num_local_cells = MIN(num_cells_per_process,
160 *num_cells - local_start);
161
162 // mask for required vertices and cells
163 int * required_vertices = xcalloc(*num_vertices, sizeof(*required_vertices));
164 int * required_cells = xcalloc(*num_cells, sizeof(*required_cells));
165
166 int offset = 0;
167 for (int i = 0; i < local_start; ++i)
168 offset += (*num_vertices_per_cell)[i];
169
170 // mark all local cells and their vertices as required
171 for (int i = local_start; i < local_start + num_local_cells; ++i) {
172
173 required_cells[i] = 2;
174 for (int j = 0; j < (*num_vertices_per_cell)[i]; ++j)
175 required_vertices[(*cell_to_vertex)[offset+j]] = 2;
176
177 offset += (*num_vertices_per_cell)[i];
178 }
179
180 // mark all halo cells as required
181 offset = 0;
182 for (int i = 0; i < *num_cells; ++i) {
183
184 if (!required_cells[i]) {
185
186 for (int j = 0; j < (*num_vertices_per_cell)[i]; ++j) {
187
188 if (required_vertices[(*cell_to_vertex)[offset+j]]) {
189
190 required_cells[i] = 1;
191 break;
192 }
193 }
194 }
195 offset += (*num_vertices_per_cell)[i];
196 }
197
198 // mark all halo vertices as required
199 offset = 0;
200 for (int i = 0; i < *num_cells; ++i) {
201
202 if (required_cells[i] == 1) {
203
204 for (int j = 0; j < (*num_vertices_per_cell)[i]; ++j) {
205
206 if (!required_vertices[(*cell_to_vertex)[offset+j]]) {
207
208 required_vertices[(*cell_to_vertex)[offset+j]] = 1;
209 }
210 }
211 }
212
213 offset += (*num_vertices_per_cell)[i];
214 }
215
216 // count the number of cells and vertices
217 int part_num_vertices = 0;
218 int part_num_cells = 0;
219 for (int i = 0; i < *num_vertices; ++i)
220 if (required_vertices[i])
221 part_num_vertices++;
222 for (int i = 0; i < *num_cells; ++i)
223 if(required_cells[i])
224 part_num_cells++;
225
226 *global_cell_id = xmalloc(part_num_cells * sizeof(**global_cell_id));
227 *cell_core_mask = xmalloc(part_num_cells * sizeof(**cell_core_mask));
228 *global_corner_id = xmalloc(part_num_vertices * sizeof(**global_corner_id));
229 *corner_core_mask = xmalloc(part_num_vertices * sizeof(**corner_core_mask));
230
231 // generate final vertex data
232 part_num_vertices = 0;
233 int * global_to_local_vertex = xmalloc(*num_vertices * sizeof(*global_to_local_vertex));
234 for (int i = 0; i < *num_vertices; ++i) {
235
236 if (required_vertices[i]) {
237
238 (*global_corner_id)[part_num_vertices] = i;
239 (*corner_core_mask)[part_num_vertices] = required_vertices[i] == 2;
240 (*x_vertices)[part_num_vertices] = (*x_vertices)[i];
241 (*y_vertices)[part_num_vertices] = (*y_vertices)[i];
242 global_to_local_vertex[i] = part_num_vertices;
243 part_num_vertices++;
244 }
245 }
246
247 *x_vertices = xrealloc(*x_vertices, part_num_vertices * sizeof(**x_vertices));
248 *y_vertices = xrealloc(*y_vertices, part_num_vertices * sizeof(**y_vertices));
249 *num_vertices = part_num_vertices;
250 free(required_vertices);
251
252 // generate final cell data
253 int num_cell_vertex_dependencies = 0;
254 part_num_cells = 0;
255 offset = 0;
256 for (int i = 0; i < *num_cells; ++i) {
257
258 if (required_cells[i]) {
259
260 (*global_cell_id)[part_num_cells] = i;
261 (*cell_core_mask)[part_num_cells] = required_cells[i] == 2;
262 (*x_cells)[part_num_cells] = (*x_cells)[i];
263 (*y_cells)[part_num_cells] = (*y_cells)[i];
264 (*cell_mask)[part_num_cells] = (*cell_mask)[i];
265
266 for (int j = 0; j < (*num_vertices_per_cell)[i]; ++j)
267 (*cell_to_vertex)[num_cell_vertex_dependencies++] =
268 global_to_local_vertex[(*cell_to_vertex)[offset+j]];
269
270 (*num_vertices_per_cell)[part_num_cells] = (*num_vertices_per_cell)[i];
271
272 part_num_cells++;
273 }
274
275 offset += (*num_vertices_per_cell)[i];
276 }
277
278 *x_cells = xrealloc(*x_cells, part_num_cells * sizeof(**x_cells));
279 *y_cells = xrealloc(*y_cells, part_num_cells * sizeof(**y_cells));
280 *cell_mask = xrealloc(*cell_mask, part_num_cells * sizeof(**cell_mask));
281
282 *num_vertices_per_cell = xrealloc(*num_vertices_per_cell, part_num_cells *
283 sizeof(**num_vertices_per_cell));
284 *cell_to_vertex = xrealloc(*cell_to_vertex, num_cell_vertex_dependencies *
285 sizeof(**cell_to_vertex));
286 *num_cells = part_num_cells;
287 free(required_cells);
288 free(global_to_local_vertex);
289}
290
291
293 char const * filename) {
294
295 int nbr_vertices;
296 int nbr_cells;
297 int * num_vertices_per_cell = NULL;
298 int * cell_to_vertex = NULL;
299
300 int * cell_mask = NULL;
301
302 double * x_vertices = NULL;
303 double * y_vertices = NULL;
304 double * x_cells = NULL;
305 double * y_cells = NULL;
306
307 yac_read_mpiom_grid_information(filename, &nbr_vertices, &nbr_cells,
309 &x_vertices, &y_vertices,
310 &x_cells, &y_cells, &cell_mask);
311
312 struct yac_basic_grid_data grid =
314 (size_t)nbr_vertices, (size_t)nbr_cells, num_vertices_per_cell,
315 x_vertices, y_vertices, cell_to_vertex);
317 free(x_vertices);
318 free(y_vertices);
319 free(cell_to_vertex);
320 free(x_cells);
321 free(y_cells);
322 free(cell_mask);
323
324 return grid;
325}
326
328 char const * filename, char const * gridname) {
329
330 return
332 gridname, yac_read_mpiom_basic_grid_data(filename));
333}
334
335/* ---------------------------------------------------------------- */
336
338 int ncid, double ** vertex_lon, double ** vertex_lat,
339 size_t * nbr_vertices) {
340
341 int glon_id;
342 int glat_id;
343
344 // get variable ids
345 if (nc_inq_varid(ncid, "lon_bounds", &glon_id) != NC_NOERR)
346 yac_nc_inq_varid(ncid, "lon_bnds", &glon_id);
347
348 if (nc_inq_varid(ncid, "lat_bounds", &glat_id) != NC_NOERR)
349 yac_nc_inq_varid(ncid, "lat_bnds", &glat_id);
350
351 int glon_ndims;
352 int glat_ndims;
353
354 // get number of dimensions
355 YAC_HANDLE_ERROR(nc_inq_varndims(ncid, glon_id, &glon_ndims));
356 YAC_HANDLE_ERROR(nc_inq_varndims(ncid, glat_id, &glat_ndims));
358 (glon_ndims == 3) && (glat_ndims == 3),
359 "ERROR(get_mpiom_vertices): "
360 "number of dimensions for lon and lat does not match")
361
362 int glon_dimids[3];
363 int glat_dimids[3];
364
365 // get dimension ids
366 YAC_HANDLE_ERROR(nc_inq_vardimid(ncid, glon_id, glon_dimids));
367 YAC_HANDLE_ERROR(nc_inq_vardimid(ncid, glat_id, glat_dimids));
368
369 size_t dimlen_lat[3];
370 size_t dimlen_lon[3];
371
372 // get size of arrays
373 for (int i = 0; i < 3; ++i) {
374 YAC_HANDLE_ERROR(nc_inq_dimlen(ncid, glon_dimids[i], &dimlen_lon[i]));
375 YAC_HANDLE_ERROR(nc_inq_dimlen(ncid, glat_dimids[i], &dimlen_lat[i]));
377 dimlen_lon[i] == dimlen_lat[i],
378 "ERROR(get_mpiom_vertices): dimlen_lon[i] != dimlen_lat[i]")
379 }
381 dimlen_lon[2] == 4, "ERROR(get_mpiom_vertices): dimlen_lon[2] has to be 4")
382
383 *nbr_vertices = dimlen_lon[0] * dimlen_lon[1] * dimlen_lon[2];
384
385 *vertex_lon = xmalloc(*nbr_vertices * sizeof(**vertex_lon));
386 *vertex_lat = xmalloc(*nbr_vertices * sizeof(**vertex_lat));
387
388 YAC_HANDLE_ERROR(nc_get_var_double(ncid, glon_id, *vertex_lon));
389 YAC_HANDLE_ERROR(nc_get_var_double(ncid, glat_id, *vertex_lat));
390}
391
392/* ---------------------------------------------------------------- */
393
394static void get_mpiom_cell_center ( int ncid,
395 double **cell_lon,
396 double **cell_lat,
397 size_t *nbr_cells ) {
398
399 int glon_id;
400 int glat_id;
401
402 // get variable ids
403 yac_nc_inq_varid (ncid, "lon", &glon_id);
404 yac_nc_inq_varid (ncid, "lat", &glat_id);
405
406 int glon_ndims;
407 int glat_ndims;
408
409 // get number of dimensions
410 YAC_HANDLE_ERROR(nc_inq_varndims(ncid, glon_id, &glon_ndims));
411 YAC_HANDLE_ERROR(nc_inq_varndims(ncid, glat_id, &glat_ndims));
413 (glon_ndims == 2) && (glat_ndims == 2),
414 "ERROR(get_mpiom_cell_center): (glon_ndims != 2) || (glat_ndims != 2)")
415
416 int glon_dimids[2];
417 int glat_dimids[2];
418
419 // get dimension ids
420 YAC_HANDLE_ERROR(nc_inq_vardimid(ncid, glon_id, glon_dimids));
421 YAC_HANDLE_ERROR(nc_inq_vardimid(ncid, glat_id, glat_dimids));
422
423 size_t dimlen_lon[2];
424 size_t dimlen_lat[2];
425
426 // get size of arrays
427 for (int i = 0; i < 2; ++i) {
428 YAC_HANDLE_ERROR(nc_inq_dimlen(ncid, glon_dimids[i], &dimlen_lon[i]));
429 YAC_HANDLE_ERROR(nc_inq_dimlen(ncid, glat_dimids[i], &dimlen_lat[i]));
431 dimlen_lon[i] == dimlen_lat[i],
432 "ERROR(get_mpiom_cell_center): dimlen_lon[i] != dimlen_lat[i]")
433 }
434
435 *nbr_cells = dimlen_lon[0] * dimlen_lon[1];
436
437 *cell_lon = xmalloc(*nbr_cells * sizeof(**cell_lon));
438 *cell_lat = xmalloc(*nbr_cells * sizeof(**cell_lat));
439
440 YAC_HANDLE_ERROR(nc_get_var_double(ncid, glon_id, *cell_lon));
441 YAC_HANDLE_ERROR(nc_get_var_double(ncid, glat_id, *cell_lat));
442}
443
444/* ---------------------------------------------------------------- */
445
447 int ncid, int **cell_mask, size_t *nbr_cells ) {
448
449 int status;
450
451 // get variable id
452 int mask_id;
453 status = nc_inq_varid (ncid, "var1", &mask_id);
454 if (status != NC_NOERR) yac_nc_inq_varid(ncid, "zo", &mask_id);
455
456 // get number of mask dimensions
457 int mask_ndims;
458 YAC_HANDLE_ERROR(nc_inq_varndims(ncid, mask_id, &mask_ndims));
460 (mask_ndims == 3) || (mask_ndims == 4),
461 "ERROR(get_mpiom_cell_mask): invalid number of dimensions")
462
463 // get dimension ids
464 int mask_dimids[4];
465 YAC_HANDLE_ERROR(nc_inq_vardimid(ncid, mask_id, mask_dimids));
466
467 // get size of mask dimensions
468 size_t dimlen_mask[4] = {0, 0, 0, 0};
469 for (int i = 0; i < mask_ndims; ++i)
470 YAC_HANDLE_ERROR(nc_inq_dimlen(ncid, mask_dimids[i], &dimlen_mask[i]));
471
472 *nbr_cells = dimlen_mask[mask_ndims - 2] * dimlen_mask[mask_ndims - 1];
473
474 int * temp_cell_mask = xmalloc(*nbr_cells * sizeof(*temp_cell_mask));
475
476 YAC_HANDLE_ERROR(nc_get_var_int(ncid, mask_id, temp_cell_mask));
477
478 if (mask_ndims == 4)
479 for (size_t i = 0; i < *nbr_cells; ++i)
480 temp_cell_mask[i] = !(temp_cell_mask[i] == 1);
481
482 *cell_mask = temp_cell_mask;
483}
484
485/* ---------------------------------------------------------------- */
486
487struct point_with_index {
488
489 int32_t lon, lat;
490 int64_t i;
491};
492
493static inline int compare_point_with_index(const void * a,const void * b) {
494
495 return ( (*(int64_t*)a < *(int64_t*)b) - (*(int64_t*)a > *(int64_t*)b) );
496}
497
499 double * temp_coords_lon, double * temp_coords_lat, int * temp_coords_mask,
500 size_t * temp_nbr_coords, int * old_to_new_id) {
501
502 struct point_with_index * sort_array =
503 xmalloc(*temp_nbr_coords * sizeof(*sort_array));
504
505 double const scale = (double)(2 << 20);
506 int32_t const periode = (int32_t)(360.0 * scale);
507
508 for (size_t i = 0; i < *temp_nbr_coords; ++i) {
509
510 int32_t curr_lon, curr_lat;
511
512 curr_lon = (int32_t)(temp_coords_lon[i] * scale);
513 curr_lat = (int32_t)(temp_coords_lat[i] * scale);
514
515 if (curr_lon <= 0)
516 curr_lon = curr_lon - (((curr_lon - periode) / periode) * periode);
517 else if (curr_lon > periode)
518 curr_lon = curr_lon - ((curr_lon / periode) * periode);
519
520 sort_array[i].lon = curr_lon;
521 sort_array[i].lat = curr_lat;
522
523 sort_array[i].i = i;
524 }
525
526 qsort(sort_array, *temp_nbr_coords, sizeof(*sort_array),
528
529 old_to_new_id[sort_array[0].i] = 1;
530
531 int last_unique_idx = (int)(sort_array[0].i);
532
533 for (size_t i = 1; i < *temp_nbr_coords; ++i) {
534
535 if (compare_point_with_index(sort_array + i - 1, sort_array + i)) {
536
537 old_to_new_id[sort_array[i].i] = 1;
538 last_unique_idx = (int)(sort_array[i].i);
539
540 } else {
541
542 old_to_new_id[sort_array[i].i] = -last_unique_idx;
543 }
544 }
545
546 free(sort_array);
547
548 size_t new_nbr_coords = 0;
549
550 for (size_t i = 0; i < *temp_nbr_coords; ++i) {
551
552 if (old_to_new_id[i] == 1) {
553
554 temp_coords_lon[new_nbr_coords] = temp_coords_lon[i];
555 temp_coords_lat[new_nbr_coords] = temp_coords_lat[i];
556
557 if (temp_coords_mask != NULL)
558 temp_coords_mask[new_nbr_coords] = temp_coords_mask[i];
559
560 old_to_new_id[i] = new_nbr_coords;
561
562 new_nbr_coords++;
563 }
564 }
565
566 for (size_t i = 0; i < *temp_nbr_coords; ++i)
567 if (old_to_new_id[i] <= 0)
568 old_to_new_id[i] = old_to_new_id[-old_to_new_id[i]];
569
570 *temp_nbr_coords = new_nbr_coords;
571}
572
573#else
574
576 char const * filename) {
577
578 UNUSED(filename);
579 die(
580 "ERROR(yac_read_mpiom_basic_grid_data): "
581 "YAC is built without the NetCDF support");
582
583 return
585 (size_t[]){0,0}, (int[]){0,0}, NULL, NULL);
586}
587
589 char const * filename, char const * gridname) {
590
591 UNUSED(filename);
592 UNUSED(gridname);
593 die(
594 "ERROR(yac_read_mpiom_basic_grid): "
595 "YAC is built without the NetCDF support");
596
597 return NULL;
598}
599
600void yac_read_mpiom_grid_information(const char * filename, int * nbr_vertices,
601 int * nbr_cells, int ** num_vertices_per_cell,
602 int ** cell_to_vertex, double ** x_vertices,
603 double ** y_vertices, double ** x_cells,
604 double ** y_cells, int ** cell_mask) {
605
606 UNUSED(filename);
607 UNUSED(nbr_vertices);
608 UNUSED(nbr_cells);
609 UNUSED(num_vertices_per_cell);
610 UNUSED(cell_to_vertex);
611 UNUSED(x_vertices);
612 UNUSED(y_vertices);
613 UNUSED(x_cells);
614 UNUSED(y_cells);
615 UNUSED(cell_mask);
616 die(
617 "ERROR(yac_read_mpiom_grid_information): "
618 "YAC is built without the NetCDF support");
619}
620
622 const char * filename, int * num_vertices, int * num_cells,
623 int ** num_vertices_per_cell, int ** cell_to_vertex, double ** x_vertices,
624 double ** y_vertices, double ** x_cells,
625 double ** y_cells, int ** global_cell_id,
626 int ** cell_mask, int ** cell_core_mask,
627 int ** global_corner_id, int ** corner_core_mask, int rank, int size) {
628
629 UNUSED(filename);
630 UNUSED(num_vertices);
631 UNUSED(num_cells);
632 UNUSED(num_vertices_per_cell);
633 UNUSED(cell_to_vertex);
634 UNUSED(x_vertices);
635 UNUSED(y_vertices);
636 UNUSED(x_cells);
637 UNUSED(y_cells);
638 UNUSED(global_cell_id);
639 UNUSED(cell_mask);
640 UNUSED(cell_core_mask);
641 UNUSED(global_corner_id);
642 UNUSED(corner_core_mask);
643 UNUSED(rank);
644 UNUSED(size);
645 die(
646 "ERROR(yac_read_part_mpiom_grid_information): "
647 "YAC is built without the NetCDF support");
648}
649
650#endif // YAC_NETCDF_ENABLED
struct yac_basic_grid * yac_basic_grid_new(char const *name, struct yac_basic_grid_data grid_data)
Definition basic_grid.c:50
struct yac_basic_grid_data yac_generate_basic_grid_data_reg_2d(size_t nbr_vertices[2], int cyclic[2], double *lon_vertices, double *lat_vertices)
Definition grid_reg2d.c:65
struct yac_basic_grid_data yac_generate_basic_grid_data_unstruct(size_t nbr_vertices, size_t nbr_cells, int *num_vertices_per_cell, double *x_vertices, double *y_vertices, int *cell_to_vertex)
#define UNUSED(x)
Definition core.h:73
#define YAC_RAD
Definition geometry.h:30
void yac_nc_inq_varid(int ncid, char const *name, int *varidp)
Definition io_utils.c:411
void yac_nc_open(const char *path, int omode, int *ncidp)
Definition io_utils.c:350
#define YAC_HANDLE_ERROR(exp)
Definition io_utils.h:30
#define xrealloc(ptr, size)
Definition ppm_xfuncs.h:67
#define xcalloc(nmemb, size)
Definition ppm_xfuncs.h:64
#define xmalloc(size)
Definition ppm_xfuncs.h:66
void yac_read_mpiom_grid_information(const char *filename, int *nbr_vertices, int *nbr_cells, int **num_vertices_per_cell, int **cell_to_vertex, double **x_vertices, double **y_vertices, double **x_cells, double **y_cells, int **cell_mask)
static void get_mpiom_cell_center(int ncid, double **cell_lon, double **cell_lat, size_t *nbr_cells)
static int compare_point_with_index(const void *a, const void *b)
static void get_mpiom_cell_mask(int ncid, int **cell_mask, size_t *nbr_cells)
static void get_mpiom_vertices(int ncid, double **vertex_lon, double **vertex_lat, size_t *nbr_vertices)
struct yac_basic_grid * yac_read_mpiom_basic_grid(char const *filename, char const *gridname)
static void remove_duplicated_coords(double *temp_coords_lon, double *temp_coords_lat, int *temp_coords_mask, size_t *temp_nbr_coords, int *old_to_new_id)
void yac_read_part_mpiom_grid_information(const char *filename, int *num_vertices, int *num_cells, int **num_vertices_per_cell, int **cell_to_vertex, double **x_vertices, double **y_vertices, double **x_cells, double **y_cells, int **global_cell_id, int **cell_mask, int **cell_core_mask, int **global_corner_id, int **corner_core_mask, int rank, int size)
struct yac_basic_grid_data yac_read_mpiom_basic_grid_data(char const *filename)
#define MIN(a, b)
#define die(msg)
Definition yac_assert.h:12
#define YAC_ASSERT(exp, msg)
Definition yac_assert.h:16