This is the second article in a series on Meteosat solar data from EUMETSAT. The intent is to define the basic parameters of meteorological data coverage for the State of Qatar. Specifically:
- Simple trigonometry is defined to assess the resolution of the satellite coverage area;
- A land surface analysis is conducted to visualize the geographic coordinates of the satellite pixels across the State of Qatar;
The geographic pixel coordinates of a Meteosat satellite are a series of standardized pixel locations across the Earth’s surface, created through data processing, which map to digital coordinates in the satellite data files. Clear understanding of pixel locations, both geographic and digital, is essential to solar data extractions. Each pixel contains high-frequency, long-term time-series data. Data available by pixel includes land surface radiation, temperature, barometric pressure, wind speeds, water vapor, cloud mask, and basic aerosol information. Data downloads can be secured from the Satellite Application Facility (SAF), EUMETSAT, Lisboa Portugal.
Spatial Resolution
Let:
Re = the Earth’s radius = 6378.14 km
d.lat(i, j) = change in lat for pixel i,j
d.lon(i,j) = change in lon for pixel i,j
lat = latitude in radians of pixel i,j
Spatial resolution in the North-South direction (ResNS) and the West-East direction (ResWE) are given by:
ResNS = Re * d.lat
ResWE = Re * d.lon * cos(lat)
Hence, Meteosat solar data resolution will vary based on the curvature of the Earth and given a specific longitude and latitude:
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 |
Meteo.Res <- function(d.lat, d.lon, lat){ # A function to estimate Meteosat pixel resolution # convert inputs from degrees to radians d.lat <- d.lat * pi/180 d.lon <- d.lon * pi/180 lat <- lat * pi/180 # Radius of the Earth R <- 6378.14 # Spatial resolution (km) in north-south direction res.NS.km <- R * d.lat # Spatial resolution (km) in west-east direction res.WE.km <- R * d.lon * cos(lat) # Combine results res.km <- matrix(c(res.NS.km, res.WE.km), ncol=1, dimnames = list(c("res.NS", "res.WE"), "km")) return(res.km) } |
Land Surface Analysis The SAF web site has a section for registered users with static data and tools located here, and includes data in HDF5 format with the geographic coordinates of the satellite pixels. The land surface image at the start of this article was generated in R using spatial data objects and the following code:
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 |
# Header ---- # # Name: LSA-FSA-pixels.R # Title: Meteosat pixel locations, resolution, and map # # Description: Simple code to map and interogate the geographic locations (lon/lat) # of Meteosat pixels. # # Details: Meteosat pixel data can be manually downloaded from SAF servers as HDF5 # files (one for lat and one for lon). The binary files are read into R using GDAL. # Separate lon/lat data is then combined into a SpatialPoint data object, and layered # on top of SpatialLines and SpatialPolygons for the State of Qatar. A function is # developed to assess the spatial resolution of the grid for any user input location. # # Pixel data: https://landsaf.meteo.pt/auxiliarDataFiles.jsp # Data file(s): 1. HDF5_LSASAF_MSG_LAT_NAfr_4bytesPrecision # 2. HDF5_LSASAF_MSG_LON_NAfr_4bytesPrecision # # GIS shores: http://www.ngdc.noaa.gov/mgg/shorelines/gshhs.html # Data file(s): 1. gshhs_f.b # # GIS country polygons: CIA World Data Bank II # Data file(s): 1. WDBII_border_f_L1 # # GIS admin boundaries Qatar: http://www.gadm.org # Data file(s): 1. QAT_adm1 # # References: LSA-SAF Product User Manual (DSSF) # EUMETSAT MFG/MSG User Handbooks # 0.Initialize Code ----- # Load libraries and DLLs library(sp, quietly = TRUE) library(rgdal, quietly = TRUE) library(maptools, quietly = TRUE) # Configure workspace and devices options(prompt = "R> ", width = 125, digits = 6, scipen = 3) map.proj <- CRS("+proj=longlat +ellps=WGS84") op=par() # 1.Load LSA-SAF pixels ---- # Load pixel coordinates and geographical coordinates lat.mena <-readGDAL("/Users/bjh/@Data/WxData/LSA-SAF/StaticData/HDF5_LSASAF_MSG_LAT_NAfr_4bytesPrecision") lon.mena <- readGDAL("/Users/bjh/@Data/WxData/LSA-SAF/StaticData/HDF5_LSASAF_MSG_LON_NAfr_4bytesPrecision") # 2.Create Pixel SpatialPoints ---- # Get lon/lat lon.sp <- lon.mena@data$band1/10000 lat.sp <- lat.mena@data$band1/10000 mat.sp <- cbind(lon.sp, lat.sp) row.names(mat.sp) <- 1:nrow(mat.sp) # Create MENA vs Qatar boundary box mena.bbox <- matrix(c(min(lon.sp, na.rm = TRUE), max(lon.sp, na.rm = TRUE), min(lat.sp, na.rm = TRUE), max(lat.sp, na.rm = TRUE)), ncol = 2, byrow = TRUE, dimnames = list(c("x", "y"), c("min", "max"))) qatar.bbox <- matrix(c(50.7, 51.65, 24.4, 26.4), ncol = 2, byrow = TRUE, dimnames = list(c("x", "y"), c("min", "max"))) # Create MENA vs Qatar spatial points mena.sp <- SpatialPoints(mat.sp, bbox = mena.bbox, proj4string = map.proj) qatar.sp <- SpatialPoints(mat.sp, bbox = qatar.bbox, proj4string = map.proj) # 3.Load Qatar Map Data ---- # GSHHS Shores Region Map (SpatialPolygons) qatar.shores <- Rgshhs("/Users/bjh/@Data/SPdata/KSA/gshhs/gshhs_f.b", xlim = c(50.65, 51.7), ylim = c(24.4, 26.4)) q.shoreGSHHS <- qatar.shores$SP # NOAA Borders (SpatialLinesDataFrame) WDBIIinfo <- ogrInfo("/Users/bjh/@Data/SPdata/KSA/WDBII_shp/f", "WDBII_border_f_L1") WDBII_b_L1 <- readOGR("/Users/bjh/@Data/SPdata/KSA/WDBII_shp/f", "WDBII_border_f_L1") q.borderNOAA <- pruneMap(WDBII_b_L1, xlim = c(50.65, 51.7), ylim = c(24.4, 26.4)) # GADM Amin Boundries GADMinfo <- ogrInfo("/Users/bjh/@Data/SPdata/Qatar/QAT-GADM", "QAT_adm1") q.adminGADM <- readOGR("/Users/bjh/@Data/SPdata/Qatar/QAT-GADM", "QAT_adm1") # Convert to SpatialPolygonDataFrame to SpatialLineDataFrame q.adminGADM <- as(q.adminGADM, "SpatialLinesDataFrame") # 4.Create State of Qatar Pixel Map ---- M0title <- "State of Qatar" M1title <- "Land Surface Analysis: Meteosat-10 Pixel Locations" S0title <- "MSG Data Servers, Satellite Application Facility, EUMETSAT" S1title <- "GSHHS and WDBII Databases, National Geophysical Data Center, NOAA" par(mar=c(2,2,2,2)+0.01, pin=c(3.73, 7.4)) plot(q.shoreGSHHS, xaxs = "i", yaxs = "i", axes = FALSE, bg = "lightcyan", col = "papayawhip", border = "lightskyblue3") plot(q.borderNOAA, led = 1.3, col = "grey45", add = TRUE) plot(q.adminGADM, led = 0.5, col = "grey60", add = TRUE) plot(qatar.sp, col = "plum", cex = 0.45, add = TRUE) axis(1, at = c(50.8 + 0:5 * 0.2), pos = 24.4, las = 1, cex.axis = 0.5, col = "grey20", col.axis = "grey20", labels=parse(text=sp:::degreeLabelsEW(c(50.8 + 0:5 * 0.2)))) axis(2, at = c(24.5 + 0:10 * 0.2), pos = 50.647, las = 2, cex.axis = 0.5, col = "grey20", col.axis = "grey20", labels = parse(text=sp:::degreeLabelsNS(c(24.5 + 0:10 * 0.2)))) axis(3, at = c(50.8 + 0:5 * 0.2), pos = 26.2865, las=1, cex.axis = 0.5, col = "grey20", col.axis = "grey20", labels = parse(text=sp:::degreeLabelsEW(c(50.8 + 0:5 * 0.2)))) axis(4, at = c(24.5 + 0:10 * 0.2), pos = 51.697, las = 2, cex.axis = 0.5, col = "grey20", col.axis = "grey20", labels = parse(text=sp:::degreeLabelsNS(c(24.5 + 0:10 * 0.2)))) title(main = M0title, font.main = 2, line = 4.50) title(main = M1title, cex.main = 0.85, font.main = 2, line = 3.00) title(sub = S0title, cex.sub = 0.80, font.sub = 3, line = 2.75) title(sub = S1title, cex.sub = 0.80, font.sub = 3, line = 4.00) box(col = "grey20", lwd = 1.3) text(x = 51.1, y = 24.9, label = "Jarayan Al Batnah", cex = 0.6) text(x = 51.46, y = 25.0, label = "Al Wakrah", cex = 0.6) text(x = 50.95, y = 25.4, label = "Al Jumayliyah", cex = 0.6) text(x = 51.15, y = 25.78, label = "Al Ghuwayriyha", cex = 0.6) text(x = 51.431, y = 25.235, label = "Ar Rayyan", cex = 0.6) text(x = 51.58, y = 25.32, label = "Ad Dawah", cex = 0.6) text(x = 51.405, y = 25.46, label = "Umm Salal", cex = 0.6) text(x = 51.39, y = 25.7, label = "Al Khawr", cex = 0.6) text(x = 51.21, y = 26, label = "Ash Shamal", cex = 0.6) par(op) |