clear all, close all, clc
load /home/anacarinelara_old/matlab/azul.mat
load /home/anacarinelara_old/matlab/brasil.txt
addpath /home2/anacarinelara/matlab/
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
%read plon and plat from regional grid file
[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;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
profw=ncread('/work3/work6/base_de_dados/WOA18/anual/woa18_decav_s00_01.nc','depth');
profp=ncread('/work3/work6/base_de_dados/MERCATOR_PISCES/clim_PICES/netCDF/Climatologia_fosfato_PISCES_m01.nc','depth');

% -------------------------------------------------------------------------

% 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('../woa2023/temp_sig2_m01.a')
% archv_fid3=fopen('pisces/nit_sig2_m01.a')
%  for 
     INDEX=1%:length(profw)
  fseek(archv_fid2,(INDEX-1)*4*(npad+IJDM),-1);
  [var,count]=fread(archv_fid2,IJDM,'float32','ieee-be');
  x=find(var>1e10); var(x)=nan;
  woa=reshape(var,IDM,JDM)'; % em metros
  clear var
% end


tt=griddata(plon,plat,woa',lonn1,latn1);
figure,pcolor(lonn1,latn1,tt),shading interp
figure,pcolor(plon,plat,woa'),shading interp

i=near(profw,1000)
%     nn=squeeze(nn_ana(:,:,i))-squeeze(nn_lua(:,:,i));
    subplot(121)
    pcolor(plon,plat,squeeze(woa(:,:,i))), shading interp, colorbar
    hold on
    plot(brasil(:,1),brasil(:,2),'k')
    title(['NIT WOA ',num2str(profw(i)),'m'])
%     caxis([0 2000]), title('Interpolado')

i=near(profp,1000)
    subplot(122)
    pcolor(plon,plat,squeeze(pis(:,:,i))), shading interp, colorbar
        hold on
    plot(brasil(:,1),brasil(:,2),'k')
%     caxis([0 2000]), 
    title(['NIT PISCES ',num2str(profp(i)),'m'])
    
i=near(plat,-23)
%     nn=squeeze(nn_ana(:,:,i))-squeeze(nn_lua(:,:,i));
    subplot(211)
    pcolor(plon,-profw,squeeze(woa(i,:,:))'), shading interp, colorbar
    hold on
    plot(brasil(:,1),brasil(:,2),'k')
    title(['NIT WOA ',num2str(plat(i)),'°S'])
    caxis([0 35])
%     caxis([0 2000]), title('Interpolado')

    subplot(212)
    pcolor(plon,-profp,squeeze(pis(i,:,:))'), shading interp, colorbar
        hold on
    plot(brasil(:,1),brasil(:,2),'k')
%     caxis([0 2000]), 
    title(['NIT PISCES ',num2str(plat(i)),'°S']) 
    caxis([0 35])
%     xlim([min(plon) max(plon)]), ylim([min(plat) max(plat)])
%     colormap(azul), caxis([-5 5])
%     title(['Prof = ',num2str(prof(i))]);
%     pause
% end

% la1 = near(lat,min(plat)); la2 = near(lat,max(plat));
% lo1 = near(lon,min(plon)); lo2 = near(lon,max(plon)); 
% n_mean=squeeze(nanmean(nanmean(nn_orig(lo1:lo2,la1:la2,:),2),1));
% n_std=squeeze(nanstd(nanstd(nn_orig(lo1:lo2,la1:la2,:),0,2),0,1));
% 
% figure,errorbar(n_mean,-prof,n_std,'horizontal','linew',2)
% grid on
% ylim([-5500 0])
% title('Perfil médio de Fosfato')
