I am trying to read an MSG (Meteosat Second Generation) file given in the native(.nat) format using satpy. I tried using:
import numpy as np
import os
from satpy import Scene
import matplotlib.pyplot as plt
import pandas as pd
from pyproj import Proj
from pyproj import transform
file = '/data/swadhin/meteosat/l1.5_data/native/MSG1-SEVI-MSG15-0100-NA-20211201035743.145000000Z-NA.nat'
# define reader
reader = "seviri_l1b_native"
# read the file
scn = Scene(filenames = {reader:[file]})
# extract data set names
dataset_names = scn.all_dataset_names()
# print available datasets
print('n'.join(map(str, dataset_names)))
which gave me
HRV
IR_016
IR_039
IR_087
IR_097
IR_108
IR_120
IR_134
VIS006
VIS008
WV_062
WV_073
I am now trying to read the 10.8 and 12 micron bands as follows:
import matplotlib.pyplot as plt
scn.load(['IR_108'])
# Access the dataset
ir108 = scn['IR_108']
# Print the metadata
print(ir108.attrs)
which again gives me:
{'orbital_parameters': {'projection_longitude': np.float32(41.5), 'projection_latitude': 0.0, 'projection_altitude': 35785831.0, 'satellite_nominal_longitude': np.float32(41.5), 'satellite_nominal_latitude': 0.0, 'satellite_actual_longitude': 41.639685531217864, 'satellite_actual_latitude': 6.875188281971676, 'satellite_actual_altitude': 35779932.76250567}, 'units': 'K', 'wavelength': WavelengthRange(min=9.8, central=10.8, max=11.8, unit='µm'), 'standard_name': 'toa_brightness_temperature', 'platform_name': 'Meteosat-8', 'sensor': 'seviri', 'georef_offset_corrected': np.True_, 'time_parameters': {'nominal_start_time': datetime.datetime(2021, 12, 1, 3, 45), 'nominal_end_time': datetime.datetime(2021, 12, 1, 4, 0), 'observation_start_time': datetime.datetime(2021, 12, 1, 3, 45, 11, 146000), 'observation_end_time': datetime.datetime(2021, 12, 1, 3, 57, 43, 145000)}, 'start_time': datetime.datetime(2021, 12, 1, 3, 45), 'end_time': datetime.datetime(2021, 12, 1, 4, 0), 'reader': 'seviri_l1b_native', 'area': Area ID: msg_seviri_iodc_3km
Description: MSG SEVIRI Indian Ocean Data Coverage service area definition with 3 km resolution
Projection: {'a': '6378169', 'h': '35785831', 'lon_0': '41.5', 'no_defs': 'None', 'proj': 'geos', 'rf': '295.488065897014', 'type': 'crs', 'units': 'm', 'x_0': '0', 'y_0': '0'}
Number of columns: 3712
Number of rows: 3712
Area extent: (np.float32(5567248.0), np.float32(5570248.5), np.float32(-5570248.5), np.float32(-5567248.0)), 'name': 'IR_108', 'resolution': 3000.403165817, 'calibration': 'brightness_temperature', 'modifiers': (), '_satpy_id': DataID(name='IR_108', wavelength=WavelengthRange(min=9.8, central=10.8, max=11.8, unit='µm'), resolution=3000.403165817, calibration=<2>, modifiers=()), 'ancillary_variables': []}
But when I try to get the brightness temperature for a particular lat and long, I tried constructing the grid as:
# Get the brightness temperature data
brightness_temp = ir108.data.compute()
# Get the area definition (contains geolocation information)
area_def = ir108.attrs['area']
# Projection parameters from ir108.attrs
proj_params = area_def.proj_dict
# Define the projection using pyproj
proj = Proj(proj='geos', h=35785831, lon_0=41.5, a=6378169, rf=295.488065897014, units='m')
# Create a grid of pixel coordinates
x = np.linspace(area_def.area_extent[0], area_def.area_extent[2], area_def.shape[1])
y = np.linspace(area_def.area_extent[1], area_def.area_extent[3], area_def.shape[0])
x_coords, y_coords = np.meshgrid(x, y)
# Transform pixel coordinates to geographic coordinates
lons, lats = proj(x_coords, y_coords, inverse=True)
# Flatten the arrays for DataFrame creation
flat_lons = lons.flatten()
flat_lats = lats.flatten()
flat_brightness_temp = brightness_temp.flatten()
# Get the time information from the metadata
start_time = ir108.attrs['start_time']
end_time = ir108.attrs['end_time']
# Create a DataFrame
df = pd.DataFrame({
'latitude': flat_lats,
'longitude': flat_lons,
'brightness_temperature': flat_brightness_temp,
'start_time': start_time,
'end_time': end_time
})
# Display the DataFrame
print(df.head())
which gave me
latitude longitude brightness_temperature start_time
0 0.0 41.500000 NaN 2021-12-01 03:45:00
1 0.0 41.500009 NaN 2021-12-01 03:45:00
2 0.0 41.500018 NaN 2021-12-01 03:45:00
3 0.0 41.500027 NaN 2021-12-01 03:45:00
4 0.0 41.500036 NaN 2021-12-01 03:45:00
end_time
0 2021-12-01 04:00:00
1 2021-12-01 04:00:00
2 2021-12-01 04:00:00
3 2021-12-01 04:00:00
4 2021-12-01 04:00:00
Where am I doing it wrong? I am doing it for the first time and can’t figure it out.
I am uploading the file here. File link – .nat file