GRIB-213: grib_get -l does not respect land-sea mask when returning nearest grid point

This commit is contained in:
Shahram Najm 2015-12-30 18:42:14 +00:00
parent c268cc8122
commit f6368eb56e
2 changed files with 23 additions and 6 deletions

View File

@ -12,7 +12,6 @@
tempLog=temp.ls.log
rm -f $tempLog
workdir=`pwd`
cd ${data_dir}
infile=regular_gaussian_model_level.grib1
@ -79,5 +78,8 @@ ${tools_dir}grib_ls -p uuidOfVGrid test_uuid.grib2 > /dev/null
type=`${tools_dir}grib_get -wcount=1 -p typeOfLevel test_uuid.grib2`
[ "$type" = "generalVertical" ]
cd $workdir
# GRIB-213 nearest with land-sea mask
temp_ls=test.grib-213.temp
${tools_dir}grib_ls -l 85,13,1,reduced_gaussian_lsm.grib1 reduced_gaussian_surface.grib1 >$temp_ls
grep -q 'Point chosen #3 index=21 .* distance=11\.' $temp_ls
rm -f $temp_ls

View File

@ -125,6 +125,9 @@ int grib_tool_init(grib_runtime_options* options)
if (options->latlon && options->latlon_mask) {
FILE* f=NULL;
grib_handle* hh;
int idx=0, land_found=0;
double min_overall = 0.0;
int idx_overall = -1;
f=fopen(options->latlon_mask,"r");
if(!f) {
perror(options->latlon_mask);
@ -146,11 +149,23 @@ int grib_tool_init(grib_runtime_options* options)
for (i=0;i<4;i++)
if (max<options->distances[i]) {max=options->distances[i];}
min=max;
min_overall=max;
/* See GRIB-213 */
for (i=0;i<4;i++) {
if ((min >= options->distances[i]) && (options->mask_values[i] >= 0.5)) {
options->latlon_idx=i;
min = options->distances[i];
if (min_overall >= options->distances[i]) { /* find overall min and index ignoring mask */
min_overall = options->distances[i];
idx_overall = i;
}
if ((min >= options->distances[i]) && (options->mask_values[i] >= 0.5)) {
land_found=1; /* found at least one point which is land */
min = options->distances[i];
idx = i;
}
}
if (land_found) {
options->latlon_idx=idx;
} else {
options->latlon_idx=idx_overall; /* all points were sea, so pick the min overall */
}
if (options->latlon_idx<0){