Skip to content

Commit

Permalink
interp and j_start_end functions
Browse files Browse the repository at this point in the history
  • Loading branch information
mlee03 authored and mlee03 committed Nov 7, 2023
1 parent 3742def commit bd19696
Show file tree
Hide file tree
Showing 3 changed files with 116 additions and 56 deletions.
64 changes: 10 additions & 54 deletions tools/fregrid/conserve_interp.c
Original file line number Diff line number Diff line change
Expand Up @@ -137,7 +137,9 @@ void setup_conserve_interp(int ntiles_in, const Grid_config *grid_in, int ntiles
}
else {

y_min = minval_double((nx_out+1)*(ny_out+1), grid_out[n].latc);
get_jstart_jend( nx_out, ny_out, nx_in, ny_in,
grid_out[n].latc, grid_in[m].latc, &jstart, &jend, &ny_now);
/*y_min = minval_double((nx_out+1)*(ny_out+1), grid_out[n].latc);
y_max = maxval_double((nx_out+1)*(ny_out+1), grid_out[n].latc);
jstart = ny_in; jend = -1;
for(j=0; j<=ny_in; j++) for(i=0; i<=nx_in; i++) {
Expand All @@ -152,7 +154,7 @@ void setup_conserve_interp(int ntiles_in, const Grid_config *grid_in, int ntiles
}
jstart = max(0, jstart-1);
jend = min(ny_in-1, jend+1);
ny_now = jend-jstart+1;
ny_now = jend-jstart+1;*/

if(opcode & CONSERVE_ORDER1) {
mxxgrid=get_maxxgrid();
Expand Down Expand Up @@ -204,64 +206,18 @@ void setup_conserve_interp(int ntiles_in, const Grid_config *grid_in, int ntiles
} // opcode GREAT_CIRCLE or CONSERVE_ORDERs


if(nxgrid > 0) {
nxgrid_prev = interp[n].nxgrid;
interp[n].nxgrid += nxgrid;
if(nxgrid_prev == 0 ) {
interp[n].i_in = (int *)malloc(interp[n].nxgrid*sizeof(int ));
interp[n].j_in = (int *)malloc(interp[n].nxgrid*sizeof(int ));
interp[n].i_out = (int *)malloc(interp[n].nxgrid*sizeof(int ));
interp[n].j_out = (int *)malloc(interp[n].nxgrid*sizeof(int ));
interp[n].area = (double *)malloc(interp[n].nxgrid*sizeof(double));
interp[n].t_in = (int *)malloc(interp[n].nxgrid*sizeof(int ));
for(i=0; i<interp[n].nxgrid; i++) {
interp[n].t_in [i] = m;
interp[n].i_in [i] = i_in [i];
interp[n].j_in [i] = j_in [i];
interp[n].i_out[i] = i_out[i];
interp[n].j_out[i] = j_out[i];
interp[n].area[i] = xgrid_area[i];
}
if(opcode & CONSERVE_ORDER2) {
interp[n].di_in = (double *)malloc(interp[n].nxgrid*sizeof(double));
interp[n].dj_in = (double *)malloc(interp[n].nxgrid*sizeof(double));
for(i=0; i<interp[n].nxgrid; i++) {
interp[n].di_in [i] = xgrid_clon[i]/xgrid_area[i];
interp[n].dj_in [i] = xgrid_clat[i]/xgrid_area[i];
}
}
}
else {
interp[n].i_in = (int *)realloc(interp[n].i_in, interp[n].nxgrid*sizeof(int ));
interp[n].j_in = (int *)realloc(interp[n].j_in, interp[n].nxgrid*sizeof(int ));
interp[n].i_out = (int *)realloc(interp[n].i_out, interp[n].nxgrid*sizeof(int ));
interp[n].j_out = (int *)realloc(interp[n].j_out, interp[n].nxgrid*sizeof(int ));
interp[n].area = (double *)realloc(interp[n].area, interp[n].nxgrid*sizeof(double));
interp[n].t_in = (int *)realloc(interp[n].t_in, interp[n].nxgrid*sizeof(int ));
for(i=0; i<nxgrid; i++) {
interp[n].t_in [nxgrid_prev+i] = m;
interp[n].i_in [nxgrid_prev+i] = i_in [i];
interp[n].j_in [nxgrid_prev+i] = j_in [i];
interp[n].i_out[nxgrid_prev+i] = i_out[i];
interp[n].j_out[nxgrid_prev+i] = j_out[i];
interp[n].area [nxgrid_prev+i] = xgrid_area[i];
}
if(opcode & CONSERVE_ORDER2) {
interp[n].di_in = (double *)realloc(interp[n].di_in, interp[n].nxgrid*sizeof(double));
interp[n].dj_in = (double *)realloc(interp[n].dj_in, interp[n].nxgrid*sizeof(double));
for(i=0; i<nxgrid; i++) {
interp[n].di_in [i+nxgrid_prev] = xgrid_clon[i]/xgrid_area[i];
interp[n].dj_in [i+nxgrid_prev] = xgrid_clat[i]/xgrid_area[i];
}
}
}
} /* if(nxgrid>0) */
get_interp( opcode, nxgrid, interp, m, n, i_in, j_in, i_out, j_out,
xgrid_clon, xgrid_clat, xgrid_area );

malloc_xgrid_arrays(zero, &i_in, &j_in, &i_out, &j_out, &xgrid_area, &xgrid_clon , &xgrid_clat);
#pragma acc exit data delete(grid_in[m].latc, grid_in[m].lonc)

} // ntiles_in

malloc_minmaxavg_lists(zero, &out_minmaxavg_lists);
#pragma acc exit data delete(out_minmaxavg_lists)
#pragma acc exit data delete(grid_out[n].latc, grid_out[n].lonc)

} // ntimes_out

if(DEBUG) print_time("time_nxgrid", time_nxgrid);
Expand Down
97 changes: 96 additions & 1 deletion tools/fregrid/conserve_interp_util.c
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
/***********************************************************************
/***********************************************************************
* GNU Lesser General Public License
*
* This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools).
Expand Down Expand Up @@ -26,6 +26,7 @@
#include "mpp.h"
#include "mpp_io.h"
#include "read_mosaic.h"
#include "mosaic_util.h"
#include "conserve_interp_util.h"

/*******************************************************************************
Expand Down Expand Up @@ -196,3 +197,97 @@ void get_CellStruct(const int tile_in, const int nx_in, const int nxgrid, int *i
}

}
/*******************************************************************************
void get_interp
********************************************************************************/
void get_interp( const int opcode, const int nxgrid, Interp_config *interp, const int m, const int n,
const int *i_in, const int *j_in, const int *i_out, const int *j_out,
const double *xgrid_clon, const double *xgrid_clat, const double *xgrid_area )
{

int nxgrid_prev;
int i;

if(nxgrid > 0) {
nxgrid_prev = interp[n].nxgrid;
interp[n].nxgrid += nxgrid;
if(nxgrid_prev == 0 ) {
interp[n].i_in = (int *)malloc(interp[n].nxgrid*sizeof(int ));
interp[n].j_in = (int *)malloc(interp[n].nxgrid*sizeof(int ));
interp[n].i_out = (int *)malloc(interp[n].nxgrid*sizeof(int ));
interp[n].j_out = (int *)malloc(interp[n].nxgrid*sizeof(int ));
interp[n].area = (double *)malloc(interp[n].nxgrid*sizeof(double));
interp[n].t_in = (int *)malloc(interp[n].nxgrid*sizeof(int ));
for(i=0; i<interp[n].nxgrid; i++) {
interp[n].t_in [i] = m;
interp[n].i_in [i] = i_in [i];
interp[n].j_in [i] = j_in [i];
interp[n].i_out[i] = i_out[i];
interp[n].j_out[i] = j_out[i];
interp[n].area[i] = xgrid_area[i];
}
if(opcode & CONSERVE_ORDER2) {
interp[n].di_in = (double *)malloc(interp[n].nxgrid*sizeof(double));
interp[n].dj_in = (double *)malloc(interp[n].nxgrid*sizeof(double));
for(i=0; i<interp[n].nxgrid; i++) {
interp[n].di_in [i] = xgrid_clon[i]/xgrid_area[i];
interp[n].dj_in [i] = xgrid_clat[i]/xgrid_area[i];
}
}
}
else {
interp[n].i_in = (int *)realloc(interp[n].i_in, interp[n].nxgrid*sizeof(int ));
interp[n].j_in = (int *)realloc(interp[n].j_in, interp[n].nxgrid*sizeof(int ));
interp[n].i_out = (int *)realloc(interp[n].i_out, interp[n].nxgrid*sizeof(int ));
interp[n].j_out = (int *)realloc(interp[n].j_out, interp[n].nxgrid*sizeof(int ));
interp[n].area = (double *)realloc(interp[n].area, interp[n].nxgrid*sizeof(double));
interp[n].t_in = (int *)realloc(interp[n].t_in, interp[n].nxgrid*sizeof(int ));
for(i=0; i<nxgrid; i++) {
interp[n].t_in [nxgrid_prev+i] = m;
interp[n].i_in [nxgrid_prev+i] = i_in [i];
interp[n].j_in [nxgrid_prev+i] = j_in [i];
interp[n].i_out[nxgrid_prev+i] = i_out[i];
interp[n].j_out[nxgrid_prev+i] = j_out[i];
interp[n].area [nxgrid_prev+i] = xgrid_area[i];
}
if(opcode & CONSERVE_ORDER2) {
interp[n].di_in = (double *)realloc(interp[n].di_in, interp[n].nxgrid*sizeof(double));
interp[n].dj_in = (double *)realloc(interp[n].dj_in, interp[n].nxgrid*sizeof(double));
for(i=0; i<nxgrid; i++) {
interp[n].di_in [i+nxgrid_prev] = xgrid_clon[i]/xgrid_area[i];
interp[n].dj_in [i+nxgrid_prev] = xgrid_clat[i]/xgrid_area[i];
}
}
}
} /* if(nxgrid>0) */

}
/*******************************************************************************
void get_jstart_jend
********************************************************************************/
void get_jstart_jend( const int nx_out, const int ny_out, const int nx_in, const int ny_in,
const double *lat_out, const double *lat_in,
int *jstart, int *jend, int *ny_now )
{

double y_min, y_max, yy ;
int i, j;

y_min = minval_double((nx_out+1)*(ny_out+1), lat_out);
y_max = maxval_double((nx_out+1)*(ny_out+1), lat_out);
*jstart = ny_in; *jend = -1;
for(j=0; j<=ny_in; j++) for(i=0; i<=nx_in; i++) {
yy = lat_in[j*(nx_in+1)+i];
if( yy > y_min ) {
if(j < *jstart ) *jstart = j;
}
if( yy < y_max ) {
if(j > *jend ) *jend = j;
}

}
*jstart = max(0, *jstart-1);
*jend = min(ny_in-1, *jend+1);
*ny_now = *jend-*jstart+1;

}
11 changes: 10 additions & 1 deletion tools/fregrid/conserve_interp_util.h
Original file line number Diff line number Diff line change
Expand Up @@ -27,8 +27,17 @@ void read_remap_file( int ntiles_in, int ntiles_out, Grid_config *grid_out,
void malloc_xgrid_arrays( int nsize, int **i_in, int **j_in, int **i_out, int **j_out,
double **xgrid_area, double **xgrid_clon, double **xgrid_clat );

void get_CellStruct(const tile_in, const int nx_in, const int nxgrid, int *i_in, int *j_in,
void get_CellStruct(const int tile_in, const int nx_in, const int nxgrid, int *i_in, int *j_in,
double *xgrid_area, double *xgrid_clon, double *xgrid_clat,
CellStruct *cell_in);

void get_interp( const int opcode, const int nxgrid, Interp_config *interp, const int m, const int n,
const int *i_in, const int *j_in, const int *i_out, const int *j_out,
const double *xgrid_clon, const double *xgrid_clat, const double *xgrid_area ) ;

void get_jstart_jend( const int nx_out, const int ny_out, const int nx_in, const int ny_in,
const double *lat_out, const double *lat_in,
int *jstart, int *jend, int *ny_now ) ;


#endif

0 comments on commit bd19696

Please sign in to comment.