search

harphub / meteogrid

R package for working with gridded meteorological data
https://harphub.github.io/meteogrid/
Other
0 stars 6 forks source link

Inconsistencies in point interpolation #20

Open andrew-MET opened 5 months ago

andrew-MET commented 5 months ago

There appears to be inconsistencies in which points are interpolated to using point.interp() that somehow depends on the length of the lat and lon vectors given to it. My suspicion is that it's got something to do with how close the points are to each other, and whether they are offshore - though no landmask is passed. Interestingly the points that aren't interpolated to with the longer list are mostly over Denmark and have an elevation of -9999, though point.interp() isn't given this information. The example below illustrates the issue. 

library(harp)
#> Loading required package: harpCore
#> 
#> Attaching package: 'harpCore'
#> The following object is masked from 'package:stats':
#> 
#>     filter
#> Loading required package: harpIO
#> Loading required package: harpPoint
#> Loading required package: harpVis
#> Loading required package: ggplot2
#> Loading required package: shiny
#> Loading required package: harpSpatial
library(meteogrid)
library(dplyr)
#> 
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union

t2m <- read_grid(
  "https://thredds.met.no/thredds/dodsC/meps25epsarchive/2024/06/28/meps_det_sfc_20240628T00Z.ncml",
  "T2m",
  lead_time = 0, 
  file_format = "netcdf"
)
#> Reading https://thredds.met.no/thredds/dodsC/meps25epsarchive/2024/06/28/meps_det_sfc_20240628T00Z.ncml

all_stns <- point.interp(t2m, station_list$lon, station_list$lat)
#> Warning: 12031 points are outside of the domain.

dom_ext <- DomainExtent(t2m)

filtered_list <- filter(
  station_list,
  between(lat, dom_ext$latlim[1], dom_ext$latlim[2]),
  between(lon, dom_ext$lonlim[1], dom_ext$lonlim[2])
)

filtered_stns <- point.interp(t2m, filtered_list$lon, filtered_list$lat)
#> Warning: 672 points are outside of the domain.

# More points found when using the filtered stations
length(na.omit(filtered_stns))
#> [1] 1569
length(na.omit(all_stns))
#> [1] 1386

# The stations that are not interpolated to from the full station list
not_included <- anti_join(
  filtered_list[!is.na(filtered_stns), ],
  station_list[!is.na(all_stns), ]
)
#> Joining with `by = join_by(SID, lat, lon, elev, name)`

print(not_included, n = Inf)
#> # A tibble: 183 × 5
#>         SID   lat   lon  elev name                 
#>       <int> <dbl> <dbl> <dbl> <chr>                
#>   1   99020  71.5 19        0 AMI (RESEARCH VESSEL)
#>   2   99090  66    2        0 FORMER MIKE          
#>   3 6005005  57.6 10.1  -9999 ---                  
#>   4 6005009  57.4  9.76 -9999 ---                  
#>   5 6005015  57.4 10.3  -9999 ---                  
#>   6 6005031  57.3 11.0  -9999 ---                  
#>   7 6005035  57.2 10.5  -9999 ---                  
#>   8 6005040  57.0 10.2  -9999 ---                  
#>   9 6005065  56.8  9.82 -9999 ---                  
#>  10 6005070  56.7 10.2  -9999 ---                  
#>  11 6005075  56.7  9.72 -9999 ---                  
#>  12 6005081  56.9  9.40 -9999 ---                  
#>  13 6005085  57.0  9.29 -9999 ---                  
#>  14 6005089  57.2  8.96 -9999 ---                  
#>  15 6005095  57.0  8.38 -9999 ---                  
#>  16 6005105  56.8  8.73 -9999 ---                  
#>  17 6005109  56.8  9.09 -9999 ---                  
#>  18 6005135  56.3  9.62 -9999 ---                  
#>  19 6005140  56.6 10.1  -9999 ---                  
#>  20 6005150  56.5 10.4  -9999 ---                  
#>  21 6005160  56.1 10.5  -9999 ---                  
#>  22 6005165  55.8 10.6  -9999 ---                  
#>  23 6005169  55.9 10.3  -9999 ---                  
#>  24 6005185  56.1  9.79 -9999 ---                  
#>  25 6005205  56.0  9.70 -9999 ---                  
#>  26 6005220  55.7 10.0  -9999 ---                  
#>  27 6005225  55.8  9.55 -9999 ---                  
#>  28 6005269  55.9  9.02 -9999 ---                  
#>  29 6005290  56.6  8.52 -9999 ---                  
#>  30 6005295  56.5  8.12 -9999 ---                  
#>  31 6005300  56.3  8.17 -9999 ---                  
#>  32 6005305  56.2  8.47 -9999 ---                  
#>  33 6005320  55.7  8.69 -9999 ---                  
#>  34 6005325  55.8  8.20 -9999 ---                  
#>  35 6005329  55.7  8.34 -9999 ---                  
#>  36 6005345  55.3  8.74 -9999 ---                  
#>  37 6005350  55.1  8.83 -9999 ---                  
#>  38 6005355  55.0  8.67 -9999 ---                  
#>  39 6005365  54.9  9.41 -9999 ---                  
#>  40 6005375  55.1  9.68 -9999 ---                  
#>  41 6005381  55.1  9.21 -9999 ---                  
#>  42 6005395  55.3  9.67 -9999 ---                  
#>  43 6005400  55.5  9.89 -9999 ---                  
#>  44 6005405  55.6 10.4  -9999 ---                  
#>  45 6005435  55.1 10.1  -9999 ---                  
#>  46 6005440  55.1 10.4  -9999 ---                  
#>  47 6005450  54.8 10.7  -9999 ---                  
#>  48 6005455  54.9 10.7  -9999 ---                  
#>  49 6005469  55.3 10.7  -9999 ---                  
#>  50 6005499  55.4 11.9  -9999 ---                  
#>  51 6005505  55.6 11.7  -9999 ---                  
#>  52 6005510  55.5 11.2  -9999 ---                  
#>  53 6005529  55.8 11.3  -9999 ---                  
#>  54 6005535  55.7 11.4  -9999 ---                  
#>  55 6005545  55.9 11.6  -9999 ---                  
#>  56 6005575  56.0 12.3  -9999 ---                  
#>  57 6005720  55.7 12.6  -9999 ---                  
#>  58 6005735  55.7 12.6  -9999 ---                  
#>  59 6005880  55.3 12.4  -9999 ---                  
#>  60 6005889  55.2 12.1  -9999 ---                  
#>  61 6005935  55.1 11.8  -9999 ---                  
#>  62 6005945  55.0 11.5  -9999 ---                  
#>  63 6005960  54.8 11.2  -9999 ---                  
#>  64 6005970  54.7 11.4  -9999 ---                  
#>  65 6005986  55.0 12.5  -9999 ---                  
#>  66 6005994  55.2 15.0  -9999 ---                  
#>  67 6020097  57.4 10.5  -9999 ---                  
#>  68 6020099  57.4 10.5  -9999 ---                  
#>  69 6020211  57.2  9.96 -9999 ---                  
#>  70 6020212  57.1 10.0  -9999 ---                  
#>  71 6020298  57.0 10.0  -9999 ---                  
#>  72 6020304  57.0  9.95 -9999 ---                  
#>  73 6020307  57.0  9.86 -9999 ---                  
#>  74 6020309  57.1  9.91 -9999 ---                  
#>  75 6020456  57.0  9.81 -9999 ---                  
#>  76 6020458  57.0  9.82 -9999 ---                  
#>  77 6020461  57.0  9.84 -9999 ---                  
#>  78 6021192  56.6  9.04 -9999 ---                  
#>  79 6021207  56.6  9.16 -9999 ---                  
#>  80 6021288  56.5  9.39 -9999 ---                  
#>  81 6021292  56.4  9.43 -9999 ---                  
#>  82 6022061  56.5 10.1  -9999 ---                  
#>  83 6022123  56.4 10.9  -9999 ---                  
#>  84 6022321  56.2 10.2  -9999 ---                  
#>  85 6022361  56.1 10.1  -9999 ---                  
#>  86 6022419  56.2  9.58 -9999 ---                  
#>  87 6022421  56.2  9.56 -9999 ---                  
#>  88 6022554  56.1 10.1  -9999 ---                  
#>  89 6023127  55.9  9.86 -9999 ---                  
#>  90 6023157  55.7  9.61 -9999 ---                  
#>  91 6023252  55.7  9.45 -9999 ---                  
#>  92 6023261  55.7  9.54 -9999 ---                  
#>  93 6023263  55.7  9.58 -9999 ---                  
#>  94 6023294  55.6  9.72 -9999 ---                  
#>  95 6023316  55.5  9.54 -9999 ---                  
#>  96 6023319  55.5  9.47 -9999 ---                  
#>  97 6023321  55.5  9.49 -9999 ---                  
#>  98 6023325  55.5  9.46 -9999 ---                  
#>  99 6023328  55.5  9.48 -9999 ---                  
#> 100 6023334  55.5  9.31 -9999 ---                  
#> 101 6023339  55.4  9.28 -9999 ---                  
#> 102 6024101  56.4  8.60 -9999 ---                  
#> 103 6024292  56.1  8.94 -9999 ---                  
#> 104 6025171  55.5  8.43 -9999 ---                  
#> 105 6026071  55.4  9.50 -9999 ---                  
#> 106 6026091  55.2  9.51 -9999 ---                  
#> 107 6026376  54.9  8.85 -9999 ---                  
#> 108 6026481  54.9  9.80 -9999 ---                  
#> 109 6028181  55.4 10.3  -9999 ---                  
#> 110 6028182  55.4 10.4  -9999 ---                  
#> 111 6028183  55.4 10.4  -9999 ---                  
#> 112 6028184  55.4 10.4  -9999 ---                  
#> 113 6028186  55.4 10.4  -9999 ---                  
#> 114 6028453  55.1 10.7  -9999 ---                  
#> 115 6028461  55.1 10.6  -9999 ---                  
#> 116 6028503  54.9 10.4  -9999 ---                  
#> 117 6029009  56.0 11.3  -9999 ---                  
#> 118 6029041  55.7 11.7  -9999 ---                  
#> 119 6029114  55.7 11.0  -9999 ---                  
#> 120 6029142  55.7 11.1  -9999 ---                  
#> 121 6029354  55.4 11.3  -9999 ---                  
#> 122 6029358  55.4 11.3  -9999 ---                  
#> 123 6029387  55.3 11.2  -9999 ---                  
#> 124 6030014  56.1 12.5  -9999 ---                  
#> 125 6030029  56.0 12.6  -9999 ---                  
#> 126 6030031  56.0 12.6  -9999 ---                  
#> 127 6030131  55.8 12.1  -9999 ---                  
#> 128 6030168  55.9 12.3  -9999 ---                  
#> 129 6030184  55.9 12.4  -9999 ---                  
#> 130 6030191  55.8 12.4  -9999 ---                  
#> 131 6030201  55.8 12.6  -9999 ---                  
#> 132 6030208  55.8 12.6  -9999 ---                  
#> 133 6030218  55.8 12.5  -9999 ---                  
#> 134 6030222  55.7 12.5  -9999 ---                  
#> 135 6030231  55.8 12.5  -9999 ---                  
#> 136 6030232  55.7 12.5  -9999 ---                  
#> 137 6030233  55.7 12.6  -9999 ---                  
#> 138 6030234  55.7 12.6  -9999 ---                  
#> 139 6030235  55.7 12.6  -9999 ---                  
#> 140 6030236  55.8 12.6  -9999 ---                  
#> 141 6030237  55.8 12.5  -9999 ---                  
#> 142 6030242  55.8 12.4  -9999 ---                  
#> 143 6030254  55.8 12.5  -9999 ---                  
#> 144 6030257  55.7 12.5  -9999 ---                  
#> 145 6030277  55.8 12.1  -9999 ---                  
#> 146 6030279  55.7 12.1  -9999 ---                  
#> 147 6030294  55.7 12.2  -9999 ---                  
#> 148 6030313  55.7 12.6  -9999 ---                  
#> 149 6030314  55.6 12.5  -9999 ---                  
#> 150 6030316  55.8 12.3  -9999 ---                  
#> 151 6030317  55.7 12.4  -9999 ---                  
#> 152 6030318  55.6 12.5  -9999 ---                  
#> 153 6030319  55.6 12.5  -9999 ---                  
#> 154 6030321  55.7 12.5  -9999 ---                  
#> 155 6030325  55.7 12.5  -9999 ---                  
#> 156 6030326  55.7 12.5  -9999 ---                  
#> 157 6030348  55.6 12.6  -9999 ---                  
#> 158 6030351  55.6 12.6  -9999 ---                  
#> 159 6030352  55.6 12.6  -9999 ---                  
#> 160 6030353  55.6 12.7  -9999 ---                  
#> 161 6030381  55.7 12.5  -9999 ---                  
#> 162 6030383  55.6 12.4  -9999 ---                  
#> 163 6030384  55.6 12.4  -9999 ---                  
#> 164 6030386  55.7 12.3  -9999 ---                  
#> 165 6030388  55.7 12.3  -9999 ---                  
#> 166 6030395  55.6 12.3  -9999 ---                  
#> 167 6030404  55.6 12.1  -9999 ---                  
#> 168 6030406  55.6 12.1  -9999 ---                  
#> 169 6030408  55.7 12.1  -9999 ---                  
#> 170 6030411  55.6 12.1  -9999 ---                  
#> 171 6030413  55.6 12.1  -9999 ---                  
#> 172 6030449  55.6 12.0  -9999 ---                  
#> 173 6030451  55.6 12.3  -9999 ---                  
#> 174 6030452  55.6 12.1  -9999 ---                  
#> 175 6031151  55.2 11.7  -9999 ---                  
#> 176 6031153  55.2 11.8  -9999 ---                  
#> 177 6031154  55.2 11.8  -9999 ---                  
#> 178 6031156  55.2 11.7  -9999 ---                  
#> 179 6031157  55.2 11.8  -9999 ---                  
#> 180 6031158  55.3 11.8  -9999 ---                  
#> 181 6031401  54.8 11.1  -9999 ---                  
#> 182 6031511  54.8 11.9  -9999 ---                  
#> 183 6032097  55.1 14.7  -9999 ---

# Plot the different stations - the rejected stations seem to either be in 
# Denmark or offshore. 

stations <- bind_rows(
  mutate(not_included, included = "no"),
  mutate(station_list[!is.na(all_stns), ], included = "yes")
)

ggplot(geo_reproject(stations, t2m), aes(x, y)) + 
  geom_polygon(
    aes(group = group), get_map(t2m), fill = "grey80", colour = "grey20"
  ) +
  geom_point(aes(colour = factor(included, levels = c("yes", "no")))) +
  scale_colour_manual("Included", values = c("blue", "red")) +
  coord_equal(expand = FALSE) +
  theme_harp_map()


# Zooming in on Denmark
ggplot(geo_reproject(stations, t2m), aes(x, y)) + 
  geom_polygon(
    aes(group = group), get_map(t2m), fill = "grey80", colour = "grey20"
  ) +
  geom_point(aes(colour = factor(included, levels = c("yes", "no")))) +
  scale_colour_manual("Included", values = c("blue", "red")) +
  coord_equal(xlim = c(-5e5, 0), ylim = c(-1e6, -5e5)) +
  theme_harp_map()

Created on 2024-06-28 with reprex v2.1.0