How to extract data only along a path consisting of multiple longitude and latitude pair points using MATLAB

Question

ANKAN SARKAR 2024-2-9

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2080326-how-to-extract-data-only-along-a-path-consisting-of-multiple-longitude-and-latitude-pair-points-usin

评论： ANKAN SARKAR 2024-2-15

Dear all,

I have gridded extinction coefficient data (model output) and I want to extract the data only along the CALIPSO swath so that I can validate the model extinction coefficient with the CALIPSO-retrieved extinction coefficient (Just like the attached screenshot which is taken from Kumar et al. (2014)). Now, I have read the longitude, latitude, and altitdue (pressure) from CALIPSO aerosol profile data and the size of both longitude and latitude vectors is 1*781 and size of altitude vector is 399*1. However, the model extinction coefficient is a matrix with size of 901*501*12 (lon*lat*plev). So, I want to extract the extinction coefficient value from the gridded model output (901*501*12) along only the 1*781 latitude-longitude pair points. How to do that in MATLAB? I tried the following code but it didn't work:

% Find the indices of the nearest grid points
lon_indices = arrayfun(@(lon) find(abs(longitude - lon) == min(abs(longitude - lon))), CALIPSO_lon);
lat_indices = arrayfun(@(lat) find(abs(latitude - lat) == min(abs(latitude - lat))), CALIPSO_lat);
% Extract data along the given latitude and longitude points
EXTCOFF_alongpoints = EXTCOFF(lon_indices, lat_indices, :);

where, 'longitude' and 'latitude' are the extracted longitude and latitude from model output, 'CALIPSO_lon' and 'CALIPSO_lat' are the longitudes and latitudes of CALIPSO aerosol profile data, and EXTCOFF is the model extinction coefficient matrix.

So, can anyone please guide me to solve this issue? That will be very helpful for me. Thank you for your time and consideration.

With regards,

Ankan

3 个评论
显示 1更早的评论隐藏 1更早的评论

ANKAN SARKAR 2024-2-13

Dear @Mathieu NOE,

I have attached all the required .mat files (Note: The 'EXTCOFF_model.mat' file with 901*501*12 (lon*lat*plev) was more than 5 MB in size, so I had to take mean over all the pressure levels to upload here. That's why it has only two dimension (901*501). But I want to extract from the mat file with lon*lat*plev dimensions so that I can plot vertical cross section just like the figure attached). Can you please guide me to solve this problem? That will be really helpful for me.

With regards,

Ankan

Mathieu NOE 2024-2-14

I think I need one more info : what lat and lon range do we have in the EXTCOFF_model data

请先登录，再进行评论。

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

Answer 1

Mathieu NOE 2024-2-14

在 MATLAB Online 中打开

ok, this is just to expose the principle... still I need from you the correct lon, lat range for the EXTCOFF_model

basically , the main task is realized with interp2

also note I played a bit with some data so the CALIPSO swath (line) would cross some interesting portion of the EXTCOFF_model data

also note I displayed the z data in log scale so the trace appears more clearly above the data

EXTCOFF_alongpoints data points appear as this line , which I slightly shifted (+0.15 in log10 base) towards higher Z values so the line appears better vs the background.

see the color (amplitude) along the line matches the EXTCOFF_model data at the corresponding lon , lat values

.

% load lon, lat data CALIPSO
load('lat_lon_alt_CALIPSO.mat')
%   alt          399x1                 3192  double              
%   lat            1x781               6248  double              
%   lon            1x781               6248  double      
r = 5;
CALIPSO_lon = lon(1:r:numel(lon))-3; % shifted this coordinate to make the line goes where I wanted 
CALIPSO_lat = lat(1:r:numel(lat))+0; % shifted this coordinate to make the line goes where I wanted 
% CALIPSO_alt = alt; % not used yet
clear lat lon alt
% load lon, lat data EXTCOFF (901*501*12 (lon*lat*plev) )
load('EXTCOFF_model.mat')
[m,n] = size(EXTCOFF);
% here it's the mean accross all altitudes so dim is (901*501*1 (lon*lat*plev) )
EXTCOFF_lon = linspace(63,73,m); % need to double check that !!
EXTCOFF_lat = linspace(0,40,n); % need to double check that !!
[X_EXTCOFF,Y_EXTCOFF] = meshgrid(EXTCOFF_lat,EXTCOFF_lon);
% Extract data along the given latitude and longitude points
 EXTCOFF_alongpoints = interp2(X_EXTCOFF,Y_EXTCOFF,EXTCOFF,CALIPSO_lat,CALIPSO_lon);
 % plot the data
figure(1),
imagesc(EXTCOFF_lat,EXTCOFF_lon,log10(EXTCOFF)); % NB log scale on z data for better rendering
colormap('jet');
hold on 
scatter3(CALIPSO_lat,CALIPSO_lon,EXTCOFF_alongpoints,50,log10(EXTCOFF_alongpoints)+0.15,'filled') % NB log scale on z data for better rendering
colorbar('vert');
hold off
figure(2),
plot(CALIPSO_lon,EXTCOFF_alongpoints)

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

Mathieu NOE 2024-2-15

在 MATLAB Online 中打开

hello

this is what I can now offer including the 3rd dimension , NB that the output is only a short segment because we have only 3 values of alt to use in the 3D interpolation (size of EXTCOFF array)

% load lon, lat data CALIPSO
load('lat_lon_alt_CALIPSO.mat')
%   alt          399x1                 3192  double              
%   lat            1x781               6248  double              
%   lon            1x781               6248  double     
% resample "alt" so it has the same number of elements as "lat" and "lon"
alt2 = interp1(1:numel(alt),alt,linspace(1,numel(alt),numel(lon)));
r = 8; % decimate so we go faster
CALIPSO_lon = lon(1:r:numel(lon));
CALIPSO_lat = lat(1:r:numel(lat));
CALIPSO_alt = alt2(1:r:numel(alt2)); 
CALIPSO_lon = CALIPSO_lon(:);
CALIPSO_lat = CALIPSO_lat(:);
CALIPSO_alt = CALIPSO_alt(:); % not used yet
clear lat lon alt alt2
% load lon, lat data EXTCOFF (901*501*12 (lon*lat*plev) )
load('EXTCOFF_model2.mat')
[m,n,l] = size(EXTCOFF);
load('lat_lon_alt_Model.mat')
%   Name              Size            Bytes  Class     Attributes
%   lat_Model       501x1              4008  double              
%   lon_Model       901x1              7208  double              
%   plev_Model       12x1                96  double  
EXTCOFF_lon = lon_Model; % 
EXTCOFF_lat = lat_Model; % 
EXTCOFF_alt = plev_Model(1:3); % first three pressure levels only to be consistent with 3rd dimension of EXTCOFF
[X_EXTCOFF,Y_EXTCOFF,Z_EXTCOFF] = meshgrid(EXTCOFF_lat,EXTCOFF_lon,EXTCOFF_alt);
% Extract data along the given latitude and longitude points
 EXTCOFF_alongpoints = interp3(X_EXTCOFF,Y_EXTCOFF,Z_EXTCOFF,EXTCOFF,CALIPSO_lat,CALIPSO_lon,CALIPSO_alt);
 % plot the data
figure(1),
EXTCOFF_mean = mean(EXTCOFF,3);
imagesc(EXTCOFF_lat,EXTCOFF_lon,log10(EXTCOFF_mean));
xlabel('Latitude');
ylabel('Longitude');
colormap('jet');
hold on 
scatter3(CALIPSO_lat,CALIPSO_lon,EXTCOFF_alongpoints,50,log10(EXTCOFF_alongpoints)+1,'filled')
colorbar('vert');
hold off
figure(2),
plot(CALIPSO_lon,EXTCOFF_alongpoints)

ANKAN SARKAR 2024-2-15

在 MATLAB Online 中打开

Hello @Mathieu NOE,

Thank you for your response. I did the following code for the proper allignment of data as per the latitude and longitude.

%%
load 'lat_lon_alt_CALIPSO.mat'
load 'lat_lon_alt_Model.mat'
%%
%lat, lon and alt of CALIPSO
lat1=lat';
lon1=flip(lon);
alt1=flip(alt);
[X_CALIPSO,Y_CALIPSO,Z_CALIPSO]=meshgrid(lon1,lat1,alt1);
%%
%lat, lon and alt of Model
lon_Model=lon_Model';
lat_Model=flip(lat_Model);
[X_Model,Y_Model,Z_Model]=meshgrid(lon_Model,lat_Model,plev_Model);
%%
%Model extinction coefficient
%Here, size of EXTCOFF is 901*501*12
EXTCOFF_1=EXTCOFF*10^3;
EXTCOFF_2 = permute(EXTCOFF_1, [2, 1, 3]);
EXTCOFF_3=flip(EXTCOFF_2,1);
EXTCOFF_alongpoints = interp3(X_Model,Y_Model,Z_Model,EXTCOFF_3,X_CALIPSO,Y_CALIPSO,Z_CALIPSO);

After running these lines, the sizes of X_CALIPSO,Y_CALIPSO, and Z_CALIPSO are equal i.e., 781*781*399 and sizes of X_Model,Y_Model, Z_Model and EXTCOFF_3 are equal i.e., 501*901*12 and sizes of EXTCOFF_alongpoints is 781*781*399. Is it correct now?

Now, I want to plot the vertical cross section (lat vs altitude) of model extinction coefficient just like the sample figure attached. Is it possible longitude-latitude pair ponts vs altitude vertcial cross section plot so that the cross-section comes only along the CALIPSO swath? Can you please help me to solve this problem? That will be very helpful for me. Thank you for your time and consideration.

With regards,

Ankan

请先登录，再进行评论。