clear all, close all, clc
load /home/anacarinelara_old/matlab/azul.mat
IDM=1500; %numeros de pontod idm
JDM=1100; %numeros de pontos jdm
KDM=32;  %numeros total de camadas
IJDM=IDM*JDM;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
npad=4096-mod(IJDM,4096);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%GRADE DO HYCOM%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%read plon and plat from regional grid file
grid_fid=fopen('/home2/mariela/hycom/HG_1_4_Ashwanth/topo/regional.grid.a','r'); %endereco do regional grid
[plon,count]=fread(grid_fid,IJDM,'float32','ieee-be');
fseek(grid_fid,4*(npad+IJDM),-1);
[plat,count]=fread(grid_fid,IJDM,'float32','ieee-be');
plon=reshape(plon,IDM,JDM); %plon=plon(:,1);
plat=reshape(plat,IDM,JDM); %plat=plat(1,:);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%% Lendo dx e dy
INDEX=9;  %%Para dx
fseek(grid_fid,INDEX*4*(npad+IJDM),-1);
[dx,count]=fread(grid_fid,IJDM,'float32','ieee-be');
INDEX=10;  %%Para dy
fseek(grid_fid,INDEX*4*(npad+IJDM),-1);
[dy,count]=fread(grid_fid,IJDM,'float32','ieee-be');
dx=reshape(dx,IDM,JDM)'; pm=1./dx;
dy=reshape(dy,IDM,JDM)'; pn=1./dy;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%% Lendo Batimetria      %%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
depth_fid=fopen('/home2/mariela/hycom/HG_1_4_Ashwanth/topo/regional.depth.a','r');
[batim,count]=fread(depth_fid,IJDM,'float32','ieee-be');
fseek(depth_fid,4*(npad+IJDM),-1);
y=find(batim>1e10);
batim(y)=nan;
batim=reshape(batim,IDM,JDM)';
dd=find(isnan(batim)==1);
cc=find(batim<800);  
%%%%%%%%%%%% CRIANDO MASCARA e inverso de DX e DY %%%%%%%%%%%%%
mask=zeros(JDM,IDM); mask(:,:)=1; mask(cc)=0; mask(dd)=0; 
mask(1,:)=0; mask(JDM,:)=0; mask(:,1)=0; mask(:,IDM)=0;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
archv_fid1=fopen('SCRATCH/relax_tem.a')
archv_fid2=fopen('../010_woa2023_old/SCRATCH/relax_tem.a')
% archv_fid2=fopen('../../010_woa2018_nit/SCRATCH/relax_tem.a')
% archv_fid2=fopen('woa2018_nit/temp_sig2_m01.a')
% archv_fid1=fopen('relax.0000_016_00.a')
NUM3D=6

% criando meshgrid para interpolar saida do modelo
lonn=[min(min(plon)):1/4:max(max(plon))];
latn=[min(min(plat)):1/4:max(max(plat))];
[lonn1,latn1]=meshgrid(lonn,latn); 

% archv_fid2=fopen('/home2/mariela/hycom/HG_1_4_Ashwanth/relax/relax.temp.a')

for kk=1:32
    INDEX=16; %:32

  fseek(archv_fid1,(kk-1)*4*(npad+IJDM),-1);
  [var,count]=fread(archv_fid1,IJDM,'float32','ieee-be');
  x=find(var>1e10); var(x)=nan;
  nn=reshape(var,IDM,JDM)'; % em metros
  clear var

  fseek(archv_fid2,(kk-1)*4*(npad+IJDM),-1);
  [var,count]=fread(archv_fid2,IJDM,'float32','ieee-be');
  x=find(var>1e10); var(x)=nan;
  nn2=reshape(var,IDM,JDM)'; % em metros
  clear var
  
  tt(:,:,kk)=griddata(plon,plat,nn',lonn1,latn1);
  tt2(:,:,kk)=griddata(plon,plat,nn2',lonn1,latn1);
end

for kk=22:32
figure(100),pcolor(lonn1,latn1,tt(:,:,kk)),shading interp
colorbar
title(['LOF: camada ', num2str(kk)])
% caxis([0 30])
caxis([0 5])

figure(200),pcolor(lonn1,latn1,tt2(:,:,kk)),shading interp
title(['Ashwanth: camada ', num2str(kk)])
colorbar
% caxis([0 30])
caxis([0 5])

figure(300),pcolor(lonn1,latn1,tt(:,:,kk)-tt2(:,:,kk)),shading interp
caxis([-2 2])

kk
pause
end

figure,pcolor(plon,plat,nn'),shading interp

for i=1:32
%     nn=squeeze(nn_ana(:,:,i))-squeeze(nn_lua(:,:,i));
figure(111)
    subplot(121)
    pcolor(plon,plat,squeeze(nn(:,:,i))), shading interp, colorbar
    caxis([22 30]), title('Relax Sal new')
    colorbar; colormap(jet)
    subplot(122)
%     pcolor(lon,lat,squeeze(nn_orig(:,:,i))'), shading interp, colorbar
%     caxis([0 1]), title('Original')
    pcolor(plon,plat,squeeze(ss(:,:,i))), shading interp, colorbar
    caxis([22 30]), title('Relax Sal old')
    colorbar; colormap(jet)

%     xlim([-50 -30]), ylim([-30 -5])
%     colormap(azul), caxis([-5 5])
%     title(['Prof = ',num2str(prof(i))]);

figure(333)
    pcolor(plon,plat,squeeze(ss(:,:,i))-squeeze(nn(:,:,i))), shading interp, colorbar
    caxis([-0.5 .5]), title('Diferença antigo - novo')
    colorbar; colormap(azul)
    pause
end

temp=ncread('/home2/base_de_dados/WOA2023/mensal/temperatura/woa23_decav_t01_01.nc','t_an');
lon=ncread('/home2/base_de_dados/WOA2023/mensal/temperatura/woa23_decav_t01_01.nc','lon');
lat=ncread('/home2/base_de_dados/WOA2023/mensal/temperatura/woa23_decav_t01_01.nc','lat');

figure,pcolor(lon,lat,squeeze(temp(:,:,1))')
shading interp
caxis([0 30])
