Skip to content
Snippets Groups Projects
Commit 539e8609 authored by Uwe Schulzweida's avatar Uwe Schulzweida
Browse files

Remapweight: prepare for multi remap files

parent b4e38162
No related branches found
No related tags found
1 merge request!143M214003/develop
Hide whitespace changes
Inline Side-by-side
libcdi @ af9b1e6f
Compare e890330e...af9b1e6f
Subproject commit e890330ebfcf84fbc95a6ab932302459314ca695
Subproject commit af9b1e6f4c456b9492e274c7b310661b34b02f24
src/Remapgrid.cc
+ 11
− 41
View file @ 539e8609
......@@ -119,29 +119,6 @@ store_gme_grid(Field &field, const Varray<int> &vgpm)
store_gme_grid(field.gridsize, field.vec_d, vgpm);
}
static int
set_maxremaps(int vlistID)
{
int maxRemaps = 0;
auto numZaxis = vlistNzaxis(vlistID);
for (int index = 0; index < numZaxis; ++index)
{
auto zaxisID = vlistZaxis(vlistID, index);
auto zaxisSize = zaxisInqSize(zaxisID);
if (zaxisSize > maxRemaps) maxRemaps = zaxisSize;
}
auto numVars = vlistNvars(vlistID);
if (numVars > maxRemaps) maxRemaps = numVars;
maxRemaps++;
if (Options::cdoVerbose) cdo_print("Set maxRemaps to %d", maxRemaps);
return maxRemaps;
}
template <typename T>
static void
remap_normalize_field(NormOpt normOpt, size_t gridsize, Varray<T> &array, T missval, const RemapGrid &tgtGrid)
......@@ -200,15 +177,6 @@ remap_set_fracmin(double fracMin, Field &field, RemapGrid *tgtGrid)
remap_set_fracmin(fracMin, field.gridsize, field.vec_d, field.missval, tgtGrid);
}
static void
remap_init(RemapType &remap)
{
remap.nused = 0;
remap.gridID = -1;
remap.gridsize = 0;
remap.nmiss = 0;
}
static void
remap_field(const RemapSwitches &remapSwitches, RemapType &remap, const Field &field1, Field &field2)
{
......@@ -431,7 +399,7 @@ public:
if (remapGrids[index]) vlistChangeGridIndex(vlistID2, index, gridID2);
maxRemaps = remapParams.maxRemaps;
if (maxRemaps == -1) maxRemaps = set_maxremaps(vlistID1);
if (maxRemaps == -1) maxRemaps = remap_get_max_maps(vlistID1);
if (maxRemaps < 1) cdo_abort("maxRemaps out of range (>0)!");
remaps.resize(maxRemaps);
......@@ -518,7 +486,9 @@ public:
field2.init(varList2[varID]);
if (!remapGrids[vlistGridIndex(vlistID1, var.gridID)]) { field_copy(field1, field2); }
auto gridIndex = vlistGridIndex(vlistID1, var.gridID);
if (!remapGrids[gridIndex]) { field_copy(field1, field2); }
else
{
if (mapType != RemapMethod::CONSERV_SCRIP && mapType != RemapMethod::CONSERV && var.gridType == GRID_GME)
......@@ -693,13 +663,6 @@ public:
if (doRemap && remap_genweights && remaps[0].nused == 0)
remap_print_warning(remapWeightsFile, operfunc, remaps[0].srcGrid, remaps[0].nmiss);
}
void
close()
{
cdo_stream_close(streamID2);
cdo_stream_close(streamID1);
for (int remapIndex = 0; remapIndex < numRemaps; remapIndex++)
{
......@@ -709,6 +672,13 @@ public:
remap_grid_free(remap.tgtGrid);
remap_search_free(remap.search);
}
}
void
close()
{
cdo_stream_close(streamID2);
cdo_stream_close(streamID1);
gridDestroy(gridID2);
......
src/Remapweights.cc
+ 67
− 32
View file @ 539e8609
......@@ -50,6 +50,7 @@ add_operators(void)
class ModuleRemapweights
{
RemapSwitches remapSwitches;
int numRemaps = 0;
int numNeighbors = 0;
CdoStreamID streamID1;
......@@ -57,27 +58,27 @@ class ModuleRemapweights
int gridID2;
Field field1;
bool useMask;
bool extrapolateIsSet;
bool remapExtrapolate;
bool needGradients;
int operfunc;
RemapType remap;
RemapMethod mapType;
int maxRemaps;
VarList varList1;
std::vector<bool> remapGrids;
std::vector<RemapType> remaps;
RemapParams remapParams;
RemapMethod mapType;
int remapOrder;
bool useMask;
bool extrapolateIsSet;
bool remapExtrapolate;
bool needGradients;
NormOpt normOpt;
bool remap_genweights = true;
int remapOrder;
NormOpt normOpt;
Field field1;
public:
void
......@@ -108,7 +109,10 @@ public:
if (inum < 1) cdo_abort("Number of nearest neighbors out of range (>0)!");
numNeighbors = inum;
}
else { operator_check_argc(1); }
else
{
operator_check_argc(1);
}
}
gridID2 = cdo_define_grid(cdo_operator_argv(0));
......@@ -125,6 +129,12 @@ public:
auto numRemapGrids = std::count_if(remapGrids.begin(), remapGrids.end(), [](auto flag) { return (flag == true); });
if (numRemapGrids == 0) cdo_abort("No remappable grid found!");
maxRemaps = remapParams.maxRemaps;
if (maxRemaps == -1) maxRemaps = remap_get_max_maps(vlistID1);
if (maxRemaps < 1) cdo_abort("maxRemaps out of range (>0)!");
remaps.resize(maxRemaps);
remapSwitches = remap_operfunc_to_maptype(operfunc);
if (numNeighbors) remapSwitches.numNeighbors = numNeighbors;
......@@ -159,32 +169,53 @@ public:
{
int varID, levelID;
cdo_inq_record(streamID1, &varID, &levelID);
field1.init(varList1[varID]);
auto &var = varList1[varID];
field1.init(var);
cdo_read_record(streamID1, field1);
auto nmiss1 = useMask ? field1.nmiss : 0;
auto gridID = varList1[varID].gridID;
auto missval = varList1[varID].missval;
auto gridsize = varList1[varID].gridsize;
auto gridIndex = vlistGridIndex(vlistID1, gridID);
auto gridIndex = vlistGridIndex(vlistID1, var.gridID);
if (remapGrids[gridIndex])
{
if (mapType != RemapMethod::CONSERV_SCRIP && mapType != RemapMethod::CONSERV && gridInqType(gridID) == GRID_GME)
if (mapType != RemapMethod::CONSERV_SCRIP && mapType != RemapMethod::CONSERV && var.gridType == GRID_GME)
cdo_abort("Only conservative remapping is available to remap between GME grids!");
if (gridIsCircular(gridID) && !extrapolateIsSet) remapExtrapolate = true;
if (gridIsCircular(var.gridID) && !extrapolateIsSet) remapExtrapolate = true;
remap_set_mask(var.gridsize, field1, nmiss1, var.missval, imask);
int remapIndex = -1;
for (remapIndex = numRemaps - 1; remapIndex >= 0; remapIndex--)
{
auto &remap = remaps[remapIndex];
if (var.gridID == remap.gridID)
{
if ((useMask && nmiss1 == remap.nmiss && imask == remap.srcGrid.mask) || !useMask)
{
remap.nused++;
break;
}
}
}
if (remapIndex >= 0) continue;
if (numRemaps >= maxRemaps) break;
remapIndex = numRemaps;
numRemaps++;
remap_set_mask(gridsize, field1, nmiss1, missval, imask);
auto &remap = remaps[remapIndex];
// remap.srcGrid.luse_cell_area = false;
// remap.tgtGrid.luse_cell_area = false;
auto numSearchBins = remapParams.numSearchBins;
if (gridInqType(gridID) != GRID_UNSTRUCTURED && !remapParams.numSearchBinsIsSet)
if (var.gridType != GRID_UNSTRUCTURED && !remapParams.numSearchBinsIsSet)
{
numSearchBins
= (!remapExtrapolate && mapType == RemapMethod::DISTWGT) ? 1 : remap_gen_numbins(gridInqYsize(gridID));
= (!remapExtrapolate && mapType == RemapMethod::DISTWGT) ? 1 : remap_gen_numbins(gridInqYsize(var.gridID));
}
remap_set_int(REMAP_NUM_SRCH_BINS, numSearchBins);
......@@ -192,19 +223,19 @@ public:
remap.vars.normOpt = normOpt;
if ((mapType == RemapMethod::BILINEAR || mapType == RemapMethod::BICUBIC)
&& (gridInqType(gridID) == GRID_GME || gridInqType(gridID) == GRID_UNSTRUCTURED))
&& (var.gridType == GRID_GME || var.gridType == GRID_UNSTRUCTURED))
cdo_abort("Bilinear/bicubic interpolation doesn't support unstructured source grids!");
// Initialize grid information for both grids
remap_init_grids(mapType, remapExtrapolate, gridID, remap.srcGrid, gridID2, remap.tgtGrid);
remap_init_grids(mapType, remapExtrapolate, var.gridID, remap.srcGrid, gridID2, remap.tgtGrid);
remap_search_init(mapType, remap.search, remap.srcGrid, remap.tgtGrid);
remap.gridID = gridID;
remap.gridID = var.gridID;
remap.nmiss = nmiss1;
if (gridInqType(gridID) == GRID_GME)
if (var.gridType == GRID_GME)
{
for (size_t i = 0, j = 0; i < gridsize; ++i)
for (size_t i = 0, j = 0; i < var.gridsize; ++i)
if (remap.srcGrid.vgpm[i]) imask[j++] = imask[i];
}
......@@ -230,11 +261,20 @@ public:
remap_gen_weights(remapParams.sortMode, remapSwitches, remap);
remap_write_data_scrip(cdo_get_stream_name(1), remapSwitches, remap.srcGrid, remap.tgtGrid, remap.vars);
break;
}
}
remap_write_data_scrip(cdo_get_stream_name(1), remapSwitches, remap.srcGrid, remap.tgtGrid, remap.vars);
for (int remapIndex = 0; remapIndex < numRemaps; remapIndex++)
{
auto &remap = remaps[remapIndex];
remap_vars_free(remap.vars);
remap_grid_free(remap.srcGrid);
remap_grid_free(remap.tgtGrid);
remap_search_free(remap.search);
}
}
void
......@@ -242,11 +282,6 @@ public:
{
cdo_stream_close(streamID1);
remap_vars_free(remap.vars);
remap_grid_free(remap.srcGrid);
remap_grid_free(remap.tgtGrid);
remap_search_free(remap.search);
cdo_finish();
}
};
......
src/remap_utils.cc
+ 32
− 0
View file @ 539e8609
......@@ -367,3 +367,35 @@ remap_set_grids(int vlistID, const VarList &varList)
return remapGrids;
}
int
remap_get_max_maps(int vlistID)
{
int maxRemaps = 0;
auto numZaxis = vlistNzaxis(vlistID);
for (int index = 0; index < numZaxis; ++index)
{
auto zaxisID = vlistZaxis(vlistID, index);
auto zaxisSize = zaxisInqSize(zaxisID);
if (zaxisSize > maxRemaps) maxRemaps = zaxisSize;
}
auto numVars = vlistNvars(vlistID);
if (numVars > maxRemaps) maxRemaps = numVars;
maxRemaps++;
if (Options::cdoVerbose) cdo_print("Set maxRemaps to %d", maxRemaps);
return maxRemaps;
}
void
remap_init(RemapType &remap)
{
remap.nused = 0;
remap.gridID = -1;
remap.gridsize = 0;
remap.nmiss = 0;
}
src/remap_utils.h
+ 4
− 0
View file @ 539e8609
......@@ -73,4 +73,8 @@ int remap_gen_numbins(int ysize);
std::vector<bool> remap_set_grids(int vlistID, const VarList &varList);
int remap_get_max_maps(int vlistID);
void remap_init(RemapType &remap);
#endif // REMAP_UTILS_H
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment

Imprint | Privacy Policy