How to plot satellite data?

Question

A.G. Bruno 2019-1-8

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/438737-how-to-plot-satellite-data

评论： Meg Noah 2023-4-24

I am currently working with IDL to plot spectral images of trace gases concentration from satellite observation (such as TROPOMI on S5P and IASI on Metop satellites) but I want to switch to MATLAB. In this process I am having several trouble:

I can not find a function or a tool to read native format images (.nat extension).
I read files in NetCDF or HDF5 formats but I can not find a good way to plot this data.

I attacched an example of what I am trying to obtain. If some one has an sample code I would be grateful if he shared it.

2 个评论
显示无隐藏无

Clare Davis 2019-2-1

Hello, I am having the same trouble - did you get a good answer in the end? Thanks

Meg Noah 2023-4-24

I do not have a .nat extension file. Can you post one here?

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Meg Noah 2020-1-10

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/438737-how-to-plot-satellite-data#answer_409609

编辑：Meg Noah 2020-1-10

在 MATLAB Online 中打开

Here's a solution

clc
close all
clear all
filename = 'W_XX-EUMETSAT-Darmstadt,HYPERSPECT+SOUNDING,METOPA+IASI_C_EUMP_20181108012658_62544_eps_o_l2.nc';
info = ncinfo(filename);
% to get a listing on console of variable names
varnames =  {info.Variables.Name}';
idx = find(contains(varnames,'integrated_co') & ~contains(varnames,'integrated_co2'));
% learning about the data
Attributes = info.Variables(idx).Attributes;
for i = 1:length(Attributes)
    if (isnumeric(Attributes(i).Value))
    fprintf(1,'%s = %f\n',Attributes(i).Name,Attributes(i).Value);
    else
    fprintf(1,'%s = %s\n',Attributes(i).Name,Attributes(i).Value);
    end
end
dataCO = ncread(filename,'integrated_co');
% they seem to have made the missing data value a really big value
minval = min(dataCO(:));
dataCO(dataCO == max(dataCO(:))) = -300;
maxval = max(dataCO(:));
fprintf(1,'minval = %f\n',minval);
fprintf(1,'maxval = %f\n',maxval);
minval = 2e-4;
maxval = 67e-4;
m = 63/(maxval-minval);
b = -m*minval;
scaleddataCO = uint8(m*dataCO+b);
scaleddataCO(dataCO<0) = nan; % missing data
colorvals = colormap(jet);
lat = ncread(filename,'lat');
lon = ncread(filename,'lon');
[ny,nx] = size(dataCO);
figure()
hold on
for i = 1:64
    idx = find(scaleddataCO == i);
    plot(lon(idx),lat(idx),'s', 'markerfacecolor',colorvals(i,:),'markeredgecolor','none');
end
xlim([-180 180])
ylim([-90 90]);
axis equal
axis tight
set(gca,'ytick',[-90 -60 -30 0 30 60 90]);
set(gca,'xtick',[-180 -150 -120 -90 -60 -30 0 30 60 90 120 150 180]);
grid on
box on
title('Integrated CO 10^{-4} kg/m^2');
colormap(jet)
hb = colorbar;
hb.Ticks = [0 0.2 0.4 0.6 0.8 1];
vals = linspace(2,67,6);
for j = 1:length(vals)
hb.TickLabels{j} = num2str(vals(j));
end

Here's how to plot it with map functions

hfig = figure('color','white','position',[32 200 1200 650]);
worldmap('World')
load coastlines
plotm(coastlat,coastlon,'k')
hfig.Name = 'Integrated CO';
hold on;
for i = 1:64
    idx = find(scaleddataCO == i);
    if (~isempty(idx))
    plotm(lat(idx),lon(idx),'s', 'markerfacecolor',colorvals(i,:),'markeredgecolor','none');
    end
end
title('Integrated CO 10^{-4} kg/m^2');
colormap(jet)
hb = colorbar('southoutside');
hb.Ticks = [0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1];
vals = linspace(2,67,11);
for j = 1:length(vals)
hb.TickLabels{j} = num2str(vals(j));
end

2 个评论
显示无隐藏无

Visweshwaran R 2020-1-27

@Meg Noah, the code is well written and precise. I too faced the same problem in plotting the netCDF file for soil moisture data (ASCAT). While I try to adjust this code for my datasets, I faced a few problems.

1) After printing the min and max value. you have assigned the same with 2e-4 and 67e-4.What is the need for it?

2)While calculating 'm', you used 63, why?

3)Similarly, you used for loop till 64, why that number particularly?

I know these questions may be basic but I am struck and 5 min of your time will be really helpful to me.

Thank you in advance.

Hugo 2023-4-24

What would it look like without the missing data part?

I'm having trouble adapting the code to a random satellite ncDEF file cause it doesnt plot all my data

请先登录，再进行评论。

Answer 2

KSSV 2019-1-8

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/438737-how-to-plot-satellite-data#answer_355465

Reading netCDF files is very easy in MATLAB. You need to just know ncdisp, ncread ..once you have read the data use pcolor to get the data you wanted to plot.

7 个评论
显示 5更早的评论隐藏 5更早的评论

A.G. Bruno 2019-1-8

在 MATLAB Online 中打开

clear all
close all
yyyy = '2018';
mm = '11';
dd = '08';
folder = '...'
file1 = '...'
filename = dir(strcat(folder, file1))
data = strcat(folder, filename.name)
    
[finfo outstrct]=read_nc_file_struct(data)
    
lon = outstrct.lon;
lat = outstrct.lat;
int_co = outstrct.integrated_co;
int_co = int_co';

This code is for a single file, to read .nc files I am currently try to use the routine read_nc_file_struct suggested by Sentinel staff. The output structures have the following dimensions.

outstrct = 
  struct with fields:
                          filename: '/scratch/atmcomp/shared/iasi/L2_metopa_nc/2018/11/08/W_XX-EUMETSAT-Darmstadt,HYPERSPECT+SOUNDING,METOPA+IASI_C_EUMP_20181108232058_62557_eps_o_l2.nc'
                   cloud_formation: [3×1 int16]
              pressure_levels_temp: [101×1 double]
          pressure_levels_humidity: [101×1 double]
             pressure_levels_ozone: [101×1 double]
    surface_emissivity_wavelengths: [12×1 double]
                 record_start_time: [765×1 double]
                  record_stop_time: [765×1 double]
                  degraded_ins_MDR: [765×1 double]
                 degraded_proc_MDR: [765×1 double]
                               lat: [120×765 double]
                               lon: [120×765 double]
                      solar_zenith: [120×765 double]
                  satellite_zenith: [120×765 double]
                     solar_azimuth: [120×765 double]
                 satellite_azimuth: [120×765 double]
        fg_atmospheric_temperature: [101×120×765 double]
        fg_atmospheric_water_vapor: [101×120×765 double]
              fg_atmospheric_ozone: [101×120×765 double]
            fg_surface_temperature: [120×765 double]
           atmospheric_temperature: [101×120×765 double]
           atmospheric_water_vapor: [101×120×765 double]
                 atmospheric_ozone: [101×120×765 double]
               surface_temperature: [120×765 double]
            integrated_water_vapor: [120×765 double]
                  integrated_ozone: [120×765 double]
                    integrated_n2o: [120×765 double]
                     integrated_co: [120×765 double]
                    integrated_ch4: [120×765 double]
                    integrated_co2: [120×765 double]
                surface_emissivity: [12×120×765 double]
           number_cloud_formations: [120×765 double]
            fractional_cloud_cover: [3×120×765 double]
             cloud_top_temperature: [3×120×765 double]
                cloud_top_pressure: [3×120×765 double]
                       cloud_phase: [3×120×765 double]
                  surface_pressure: [120×765 double]
                   instrument_mode: [765×1 double]
               spacecraft_altitude: [765×1 double]
                      flag_amsubad: [120×765 double]
                     flag_avhrrbad: [120×765 double]
                       flag_cdlfrm: [120×765 double]
                       flag_cldnes: [120×765 double]
                       flag_cdltst: [120×765 uint16]
                       flag_daynit: [120×765 double]
                      flag_dustcld: [120×765 double]
                      flag_fgcheck: [120×765 uint16]
                      flag_iasibad: [120×765 double]
                       flag_initia: [120×765 double]
                       flag_itconv: [120×765 double]
                      flag_landsea: [120×765 double]
                       flag_mhsbad: [120×765 double]
                        flag_numit: [120×765 double]
                       flag_nwpbad: [120×765 double]
                    flag_physcheck: [120×765 double]
                     flag_retcheck: [120×765 uint16]
                       flag_satman: [120×765 double]
                      flag_sunglnt: [120×765 double]
                       flag_thicir: [120×765 double]
                       nerr_values: [765×1 double]
                  error_data_index: [120×765 double]
                 temperature_error: [406×116×765 uint32]
                water_vapour_error: [171×116×765 uint32]
                       ozone_error: [55×116×765 uint32]
                         surface_z: [120×765 double]
                          co_qflag: [120×765 double]
                           co_bdiv: [120×765 uint32]
                           co_npca: [120×765 double]
                     co_nfitlayers: [120×765 double]
                     co_nbr_values: [765×1 double]
                         co_cp_air: [19×120×765 double]
                        co_cp_co_a: [19×120×765 double]
                           co_x_co: [19×120×765 double]
                  co_h_eigenvalues: [10×120×765 double]
                 co_h_eigenvectors: [190×120×765 double]

KSSV 2019-1-9

在 MATLAB Online 中打开

What are the dimensions of each variable?

Try:

pcolor(lon,lat,int_co') ;

I have asked you to attach the data/ file so that you can get more help straight away.

A.G. Bruno 2019-1-10

It still doesn't work.

Sorry, I did not understand that you mean the data file. I've attached a dropbox link because the file is too large to upload it.

https://www.dropbox.com/s/2dkep0bzimms66m/W_XX-EUMETSAT-Darmstadt%2CHYPERSPECT%2BSOUNDING%2CMETOPA%2BIASI_C_EUMP_20181108012658_62544_eps_o_l2.7z?dl=0

请先登录，再进行评论。