MATLAB Answers

How to set level values of a variable by comparing(subtracing) two nc files

4 ビュー (過去 30 日間)
Karthik M
Karthik M 2021 年 9 月 14 日
コメント済み: Karthik M 2021 年 9 月 14 日
Hi Freinds
I need to spot seasonal trends of chlor_a nc data for last 10 years for which I considered 2 year(2015 and 2016) nc files (Feb)month at the moment
I find error while subtracting both files to compare and identify max persistence (chlor_a conc.) hotspot of particular month during successive year . After executing need to save both files in nc format, [+ve] value means max chlor_a conc. and [-ve] value means less conc.
I have extracted the chlor_a data and have plot with respective x and y coordinates,
I used below code for both nc files
ncfile_1 = 'A20150322015059.L3m_MO_CHL_chlor_a_4km.nc';
lat_1 = ncread(ncfile_1,'lat') ;
lon_1 = ncread(ncfile_1,'lon') ;
chlor_a_1 = ncread(ncfile_1,'chlor_a');
[X_1,Y_1] = meshgrid(lon_1,lat_1) ;
nn = 1000;
xi_1 = linspace(30,100,nn) ;
yi_1 = linspace(0,30,nn) ;
[Xi_1,Yi_1] = meshgrid(xi_1,yi_1);
iwant_1 = interp2(X_1,Y_1,chlor_a_1',Xi_1,Yi_1)' ;
ind_1 = find(iwant_1>10 | iwant_1<0);
%ind = find(iwant_2 > 0 & iwant_1<10);
iwant_1(ind_1) = NaN;
%[min_val_1,min_ind_1] = min(iwant_1',[],'all','linear');
[min_val_1,min_ind_1] = min(reshape(iwant_1',[],1));
[r1_1,c1_1]=ind2sub(nn,min_ind_1);
x_min_1 = xi_1(c1_1);
y_min_1 = yi_1(r1_1);
%[max_val_1,max_ind_1] = max(iwant_1',[],'all','linear');
[max_val_1,max_ind_1] = max(reshape(iwant_1',[],1));
[r2_1,c2_1]=ind2sub(nn,max_ind_1);
x_max_1 = xi_1(c2_1);
y_max_1 = yi_1(r2_1);
figure(1)
pcolor(xi_1,yi_1,iwant_1'); shading interp;
c = colorbar;
cmap = jet(255);
cmap(:,2) = 0;
colormap(cmap)
caxis([0, 5])
hold on
plot(x_min_1,y_min_1,'k +', 'MarkerSize', 10);
plot(x_max_1,y_max_1,'m +', 'MarkerSize', 10);
hold off
legend('','min','max');
I need to get the plot as such after subtraction, please find the attached image
similary need to work on to compare weekly data
data link
Kindly help, waitng for your valuable inputs
Thanks in advance
regards
  4 件のコメント
Karthik M
Karthik M 2021 年 9 月 14 日
Dear @KSSV
without setting the index found the error image as attached. I want to set level value of variable with >10 and <0 = to ignore the unwanted values with the line in my code [ind = find(iwant_2>10 | iwant_2<0)]; as shown
While subtracting as 2 mat nc file_1 and ncfile_2 (with index limit). I am trying to obtain the sample image attached from the below code
ncfile_1 = 'A20150322015059.L3m_MO_CHL_chlor_a_4km.nc';
lat_1 = ncread(ncfile_1,'lat') ;
lon_1 = ncread(ncfile_1,'lon') ;
chlor_a_1 = ncread(ncfile_1,'chlor_a');
[X_1,Y_1] = meshgrid(lon_1,lat_1) ;
nn = 1000;
xi_1 = linspace(30,100,nn) ;
yi_1 = linspace(0,30,nn) ;
[Xi_1,Yi_1] = meshgrid(xi_1,yi_1);
iwant_1 = interp2(X_1,Y_1,chlor_a_1',Xi_1,Yi_1)' ;
ind_1 = find(iwant_1>10 | iwant_1<0);
%ind = find(iwant_2 > 0 & iwant_1<10);
iwant_1(ind_1) = NaN;
%[min_val_1,min_ind_1] = min(iwant_1',[],'all','linear');
[min_val_1,min_ind_1] = min(reshape(iwant_1',[],1));
[r1_1,c1_1]=ind2sub(nn,min_ind_1);
x_min_1 = xi_1(c1_1);
y_min_1 = yi_1(r1_1);
%[max_val_1,max_ind_1] = max(iwant_1',[],'all','linear');
[max_val_1,max_ind_1] = max(reshape(iwant_1',[],1));
[r2_1,c2_1]=ind2sub(nn,max_ind_1);
x_max_1 = xi_1(c2_1);
y_max_1 = yi_1(r2_1);
figure(1)
pcolor(xi_1,yi_1,iwant_1'); shading interp;
c = colorbar;
cmap = jet(255);
cmap(:,2) = 0;
colormap(cmap)
caxis([0, 5])
hold on
plot(x_min_1,y_min_1,'k +', 'MarkerSize', 10);
plot(x_max_1,y_max_1,'m +', 'MarkerSize', 10);
hold off
legend('','min','max');
Image = getframe(gcf);
imwrite(Image.cdata, 'feb2015.tif');
--
ncfile_2 = 'A20160322016060.L3m_MO_CHL_chlor_a_4km.nc';
lat_2 = ncread(ncfile_2,'lat') ;
lon_2 = ncread(ncfile_2,'lon') ;
chlor_a_2 = ncread(ncfile_2,'chlor_a');
[X_2,Y_2] = meshgrid(lon_2,lat_2) ;
nn = 1000;
xi_2 = linspace(30,100,nn) ;
yi_2 = linspace(0,30,nn) ;
[Xi_2,Yi_2] = meshgrid(xi_2,yi_2);
iwant_2 = interp2(X_2,Y_2,chlor_a_2',Xi_2,Yi_2)' ;
ind = find(iwant_2>10 | iwant_2<0);
%ind = find(iwant_2 > 0 & iwant_2<10);
iwant_2(ind) = NaN;
%[min_val,min_ind] = min(iwant_2',[],'all','linear');
[min_val,min_ind] = min(reshape(iwant_2',[],1));
[r1,c1]=ind2sub(nn,min_ind);
x_min = xi_2(c1);
y_min = yi_2(r1);
%[max_val,max_ind] = max(iwant_2',[],'all','linear');
[max_val,max_ind] = max(reshape(iwant_2',[],1));
[r2,c2]=ind2sub(nn,max_ind);
x_max = xi_2(c2);
y_max = yi_2(r2);
figure(2)
pcolor(xi_2,yi_2,iwant_2'); shading interp;
c = colorbar;
cmap = jet(255);
cmap(:,2) = 0;
colormap(cmap)
caxis([0, 5])
hold on
plot(x_min,y_min,'k +', 'MarkerSize', 10);
plot(x_max,y_max,'m +', 'MarkerSize', 10);
hold off
legend('','min','max');
I have tried my level best for subtracting the 2 files after indexing. I don't know how to execute
kindly help, need to execute file as sample image
thank you
regards

サインインしてコメントする。

採用された回答

KSSV
KSSV 2021 年 9 月 14 日
Your interpolation grid coordinates are completely different comapred to original grid. Read about interp2.
ncfile_1 = 'A20150322015059.L3m_MO_CHL_chlor_a_4km.nc';
lat_1 = ncread(ncfile_1,'lat') ;
lon_1 = ncread(ncfile_1,'lon') ;
chlor_a_1 = ncread(ncfile_1,'chlor_a');
[X_1,Y_1] = meshgrid(lon_1,lat_1) ;
xi_1 = linspace(min(lon_1),max(lon_1),1000) ; % Changed here
yi_1 = linspace(min(lat_1),max(lat_1),1000) ; % changed here
[Xi_1,Yi_1] = meshgrid(xi_1,yi_1);
iwant_1 = interp2(X_1,Y_1,chlor_a_1',Xi_1,Yi_1)' ;
figure(1)
pcolor(xi_1,yi_1,iwant_1'); shading interp;
c = colorbar;
cmap = jet(255);
cmap(:,2) = 0;
colormap(cmap)
caxis([0, 5])
ncfile_2 = 'A20160322016060.L3m_MO_CHL_chlor_a_4km.nc';
lat_2 = ncread(ncfile_2,'lat') ;
lon_2 = ncread(ncfile_2,'lon') ;
chlor_a_2 = ncread(ncfile_2,'chlor_a');
[X_2,Y_2] = meshgrid(lon_2,lat_2) ;
xi_2 = linspace(min(lon_2),max(lon_2),1000) ; % Changed here
yi_2 = linspace(min(lat_2),max(lat_2),1000) ; % Changed here
[Xi_2,Yi_2] = meshgrid(xi_2,yi_2);
iwant_2 = interp2(X_2,Y_2,chlor_a_2',Xi_2,Yi_2)' ;
figure(2)
pcolor(xi_2,yi_2,iwant_2'); shading interp;
c = colorbar;
cmap = jet(255);
cmap(:,2) = 0;
colormap(cmap)
caxis([0, 5])
figure(3)
diff = iwant_2-iwant_1;
%M = max(diff, [], 'all');
%M
pcolor(diff'); shading interp;
grid on;
imshow(diff);
imgdiff( 4, 3, : ) = [0 0 0];
mask = sum(imgdiff, 3 ) > 0;
h = imagesc(imgdiff);
h.AlphaData = mask;
%axis([0 30,30 100]);
c = colorbar;
cmap = jet(255);
cmap(:,2) = 0;
colormap(cmap)
caxis([0, 5])
  5 件のコメント
Karthik M
Karthik M 2021 年 9 月 14 日
Thank you @KSSV for your valuable time

サインインしてコメントする。

その他の回答 (0 件)

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by