From 8bf77373f8a547fabc20f17aee3b6d75a05b6d27 Mon Sep 17 00:00:00 2001
From: Christian Ethe <cetlod@espri-calcul-01.ipsl.fr>
Date: Mon, 10 Mar 2025 15:54:38 +0100
Subject: [PATCH] Minor improvments to create PISCES atmospheric deposition
 files for one or several childsgrid

---
 Nesting_tools/create_nesteddust.m  |  36 ++++++++-
 Nesting_tools/create_nestedndepo.m | 113 +++++++++++++++++++++++++++
 Nesting_tools/nested_clim.m        |   2 +-
 Nesting_tools/nested_dust.m        |   4 +
 Nesting_tools/nested_ndepo.m       | 118 +++++++++++++++++++++++++++++
 Nesting_tools/plot_nestndepo.m     |  83 ++++++++++++++++++++
 6 files changed, 353 insertions(+), 3 deletions(-)
 create mode 100644 Nesting_tools/create_nestedndepo.m
 create mode 100644 Nesting_tools/nested_ndepo.m
 create mode 100644 Nesting_tools/plot_nestndepo.m

diff --git a/Nesting_tools/create_nesteddust.m b/Nesting_tools/create_nesteddust.m
index ffb0ac73..3de54489 100644
--- a/Nesting_tools/create_nesteddust.m
+++ b/Nesting_tools/create_nesteddust.m
@@ -37,7 +37,7 @@ close(nc);
 Lp=L+1;
 Mp=M+1;
 nw = netcdf(dustname, 'clobber');
-%redef(nw);
+redef(nw);
 %
 %  Create dimensions
 %
@@ -70,7 +70,39 @@ nw{'dust'}.units = 'nmol Fe m-3';
 nw{'dust'}.field = ncchar('Fe Dust Deposition, scalar, series');
 nw{'dust'}.field = 'Fe Dust Deposition, scalar, series';
 
-%endef(nw);
+nw{'dustfer'} = ncdouble('dust_time', 'eta_rho', 'xi_rho');
+nw{'dustfer'}.long_name = ncchar('Fe Dust Deposition');
+nw{'dustfer'}.long_name = 'Fe Dust Deposition';
+nw{'dustfer'}.units = ncchar('Kg m-2 s-1');
+nw{'dustfer'}.units = 'Kg m-2 s-1';
+nw{'dustfer'}.field = ncchar('dustfer, scalar, series');
+nw{'dustfer'}.field = 'dustfer, scalar, series';
+%
+nw{'dustpo4'} = ncdouble('dust_time', 'eta_rho', 'xi_rho');
+nw{'dustpo4'}.long_name = ncchar('PO4 Dust Deposition');
+nw{'dustpo4'}.long_name = 'PO4 Dust Deposition';
+nw{'dustpo4'}.units = ncchar('Kg m-2 s-1');
+nw{'dustpo4'}.units = 'Kg m-2 s-1';
+nw{'dustpo4'}.field = ncchar('dustpo4, scalar, series');
+nw{'dustpo4'}.field = 'dustpo4, scalar, series';
+%
+nw{'dustsi'} = ncdouble('dust_time', 'eta_rho', 'xi_rho');
+nw{'dustsi'}.long_name = ncchar('Si Dust Deposition');
+nw{'dustsi'}.long_name = 'Si Dust Deposition';
+nw{'dustsi'}.units = ncchar('Kg m-2 s-1');
+nw{'dustsi'}.units = 'Kg m-2 s-1';
+nw{'dustsi'}.field = ncchar('dustsi, scalar, series');
+nw{'dustsi'}.field = 'dustsi, scalar, series';
+%
+nw{'solubility2'} = ncdouble('dust_time', 'eta_rho', 'xi_rho');
+nw{'solubility2'}.long_name = ncchar('Fe solubility from Mahowald');
+nw{'solubility2'}.long_name = 'Fe solubility from Mahowald';
+nw{'solubility2'}.units = ncchar('%');
+nw{'solubility2'}.units = '%';
+nw{'solubility2'}.field = ncchar('solubility2, scalar, series');
+nw{'solubility2'}.field = 'solubility2, scalar, series';
+
+endef(nw);
 
 %
 % Create global attributes
diff --git a/Nesting_tools/create_nestedndepo.m b/Nesting_tools/create_nestedndepo.m
new file mode 100644
index 00000000..18ce5517
--- /dev/null
+++ b/Nesting_tools/create_nestedndepo.m
@@ -0,0 +1,113 @@
+function  create_nestedndepo(ndeponame,parentname,grdname,title,ndepot,ndepoc)
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%
+% 	Create an empty netcdf ndepo forcing file (for PISCES)
+%       ndeponame: name of the ndepo file
+%       grdname: name of the grid file
+%       title: title in the netcdf file  
+%
+%  Further Information:  
+%  http://www.croco-ocean.org
+%  
+%  This file is part of CROCOTOOLS
+%
+%  CROCOTOOLS is free software; you can redistribute it and/or modify
+%  it under the terms of the GNU General Public License as published
+%  by the Free Software Foundation; either version 2 of the License,
+%  or (at your option) any later version.
+%
+%  CROCOTOOLS is distributed in the hope that it will be useful, but
+%  WITHOUT ANY WARRANTY; without even the implied warranty of
+%  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+%  GNU General Public License for more details.
+%
+%  You should have received a copy of the GNU General Public License
+%  along with this program; if not, write to the Free Software
+%  Foundation, Inc., 59 Temple Place, Suite 330, Boston,
+%  MA  02111-1307  USA
+%
+%  Copyright (c) 2004-2006 by Pierrick Penven 
+%  e-mail:Pierrick.Penven@ird.fr  
+%  Update : Gildas Cambon 13 Oct 2009
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+nc=netcdf(grdname);
+L=length(nc('xi_psi'));
+M=length(nc('eta_psi'));
+close(nc);
+Lp=L+1;
+Mp=M+1;
+nw = netcdf(ndeponame, 'clobber');
+redef(nw);
+%
+%  Create dimensions
+%
+nw('xi_u') = L;
+nw('eta_u') = Mp;
+nw('xi_v') = Lp;
+nw('eta_v') = M;
+nw('xi_rho') = Lp;
+nw('eta_rho') = Mp;
+nw('xi_psi') = L;
+nw('eta_psi') = M;
+nw('ndepo_time') = length(ndepot);
+%
+%  Create variables and attributes
+%
+nw{'ndepo_time'} = ncdouble('ndepo_time');
+nw{'ndepo_time'}.long_name = ncchar('ndepo time');
+nw{'ndepo_time'}.long_name = 'ndepo time';
+nw{'ndepo_time'}.units = ncchar('days');
+nw{'ndepo_time'}.units = 'days';
+nw{'ndepo_time'}.cycle_length = ndepoc;
+nw{'ndepo_time'}.field = ncchar('time, scalar, series');
+nw{'ndepo_time'}.field = 'time, scalar, series';
+
+%
+nw{'ndepo'} = ncdouble('ndepo_time', 'eta_rho', 'xi_rho');
+nw{'ndepo'}.long_name = ncchar('Nitrogen Deposition');
+nw{'ndepo'}.long_name = 'Nitrogen Deposition';
+nw{'ndepo'}.units = ncchar('KgN m-2 s-1');
+nw{'ndepo'}.units = 'KgN m-2 s-1';
+nw{'ndepo'}.field = ncchar('ndepo, scalar, series');
+nw{'ndepo'}.field = 'ndepo, scalar, series';
+%
+nw{'noyndepo'} = ncdouble('ndepo_time', 'eta_rho', 'xi_rho');
+nw{'noyndepo'}.long_name = ncchar('NOy Deposition');
+nw{'noyndepo'}.long_name = 'NOy Deposition';
+nw{'noyndepo'}.units = ncchar('KgN m-2 s-1');
+nw{'noyndepo'}.units = 'KgN m-2 s-1';
+nw{'noyndepo'}.field = ncchar('noyndepo, scalar, series');
+nw{'noyndepo'}.field = 'noyndepo, scalar, series';
+%
+nw{'nhxndepo'} = ncdouble('ndepo_time', 'eta_rho', 'xi_rho');
+nw{'nhxndepo'}.long_name = ncchar('NHx Deposition');
+nw{'nhxndepo'}.long_name = 'NHx Deposition';
+nw{'nhxndepo'}.units = ncchar('KgN m-2 s-1');
+nw{'nhxndepo'}.units = 'KgN m-2 s-1';
+nw{'nhxndepo'}.field = ncchar('nhxndepo, scalar, series');
+nw{'nhxndepo'}.field = 'nhxndepo, scalar, series';
+%
+endef(nw);
+
+%
+% Create global attributes
+%
+
+nw.title = ncchar(title);
+nw.title = title;
+nw.date = ncchar(date);
+nw.date = date;
+nw.grd_file = ncchar(grdname);
+nw.grd_file = grdname;
+nw.parent_file = ncchar(parentname);
+nw.parent_file = parentname;
+
+%
+% Write time variables
+%
+
+nw{'ndepo_time'}(:) = ndepot;
+
+
+close(nw);
+return
diff --git a/Nesting_tools/nested_clim.m b/Nesting_tools/nested_clim.m
index b3847285..a31095c2 100644
--- a/Nesting_tools/nested_clim.m
+++ b/Nesting_tools/nested_clim.m
@@ -57,7 +57,7 @@ elseif (nvar <=  41 & biol)
     disp(['Compute Biogeochemical variables type NPZD : '])
     disp(['NChlPZD or N2ChlPZD2                     '])
     disp('==========================')
-elseif (46 <= nvar <=47 & biobus)
+elseif (46 <= nvar <=47 & bioebus)
     timebiol={'no3_time';'o2_time';'chla_time';'sphyto_time';'lphyto_time';'szoo_time';'lzoo_time';'n2o_time'};
     cyclebiol={'no3_cyle';'o2_cycle';'chla_cycle';'sphyto_cycle';'lphyto_cycle';'szoo_cycle';'lzoo_cycle';'n2o_cycle'};
     namebiol={'NO3';'O2';'CHLA';'SPHYTO';'LPHYTO';'SZOO';'LZOO';'N2O'};
diff --git a/Nesting_tools/nested_dust.m b/Nesting_tools/nested_dust.m
index 98a94db7..10d34370 100644
--- a/Nesting_tools/nested_dust.m
+++ b/Nesting_tools/nested_dust.m
@@ -102,6 +102,10 @@ nc=netcdf(child_dust,'write');
 disp('dust...')
 for tindex=1:length(dustt)
   interpvar3d(np,nc,igrd_r,jgrd_r,ichildgrd_r,jchildgrd_r,'dust',mask,tindex)
+  interpvar3d(np,nc,igrd_r,jgrd_r,ichildgrd_r,jchildgrd_r,'dustfer',mask,tindex)
+  interpvar3d(np,nc,igrd_r,jgrd_r,ichildgrd_r,jchildgrd_r,'dustpo4',mask,tindex)
+  interpvar3d(np,nc,igrd_r,jgrd_r,ichildgrd_r,jchildgrd_r,'dustsi',mask,tindex)
+  interpvar3d(np,nc,igrd_r,jgrd_r,ichildgrd_r,jchildgrd_r,'solubility2',mask,tindex)
 end
 result=close(np);
 result=close(nc);
diff --git a/Nesting_tools/nested_ndepo.m b/Nesting_tools/nested_ndepo.m
new file mode 100644
index 00000000..77df9324
--- /dev/null
+++ b/Nesting_tools/nested_ndepo.m
@@ -0,0 +1,118 @@
+function nested_ndepo(child_grd,parent_ndepo,child_ndepo)
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%
+%  compute the ndepo file (PISCES biogeochemical model) 
+%  of the embedded grid
+%
+%  Further Information:  
+%  http://www.croco-ocean.org
+%  
+%  This file is part of CROCOTOOLS
+%
+%  CROCOTOOLS is free software; you can redistribute it and/or modify
+%  it under the terms of the GNU General Public License as published
+%  by the Free Software Foundation; either version 2 of the License,
+%  or (at your option) any later version.
+%
+%  CROCOTOOLS is distributed in the hope that it will be useful, but
+%  WITHOUT ANY WARRANTY; without even the implied warranty of
+%  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+%  GNU General Public License for more details.
+%
+%  You should have received a copy of the GNU General Public License
+%  along with this program; if not, write to the Free Software
+%  Foundation, Inc., 59 Temple Place, Suite 330, Boston,
+%  MA  02111-1307  USA
+%
+%  Copyright (c) 2004-2006 by Pierrick Penven 
+%  e-mail:Pierrick.Penven@ird.fr  
+%  Update : Gildas Cambon: 13 Nov 2009
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+extrapmask=0;
+%
+% Title
+%
+title=['ndepo file for the embedded grid :',child_ndepo,...
+       ' using parent ndepo file: ',parent_ndepo];
+disp(' ')
+disp(title)
+%
+% Read in the embedded grid
+%
+disp(' ')
+disp(' Read in the embedded grid...')
+nc=netcdf(child_grd);
+parent_grd=nc.parent_grid(:);
+imin=nc{'grd_pos'}(1);
+imax=nc{'grd_pos'}(2);
+jmin=nc{'grd_pos'}(3);
+jmax=nc{'grd_pos'}(4);
+refinecoeff=nc{'refine_coef'}(:);
+result=close(nc);
+nc=netcdf(parent_grd);
+Lp=length(nc('xi_rho'));
+Mp=length(nc('eta_rho'));
+if extrapmask==1
+  mask=nc{'mask_rho'}(:);
+else
+  mask=[];
+end
+result=close(nc);
+%
+% Read in the parent ndepo file
+%
+disp(' ')
+disp(' Read in the parent ndepo file...')
+nc = netcdf(parent_ndepo);
+ndepot = nc{'ndepo_time'}(:);
+ndepoc = nc{'ndepo_time'}.cycle_length(:);
+result=close(nc);
+%
+% Create the ndepo file
+%
+disp(' ')
+disp(' Create the ndepo file...')
+create_nestedndepo(child_ndepo,parent_ndepo,child_grd,title,...
+		  ndepot,ndepoc)
+%
+% parent indices
+%
+[igrd_r,jgrd_r]=meshgrid((1:1:Lp),(1:1:Mp));
+[igrd_p,jgrd_p]=meshgrid((1:1:Lp-1),(1:1:Mp-1));
+[igrd_u,jgrd_u]=meshgrid((1:1:Lp-1),(1:1:Mp));
+[igrd_v,jgrd_v]=meshgrid((1:1:Lp),(1:1:Mp-1));
+%
+% the children indices
+%
+ipchild=(imin:1/refinecoeff:imax);
+jpchild=(jmin:1/refinecoeff:jmax);
+irchild=(imin+0.5-0.5/refinecoeff:1/refinecoeff:imax+0.5+0.5/refinecoeff);
+jrchild=(jmin+0.5-0.5/refinecoeff:1/refinecoeff:jmax+0.5+0.5/refinecoeff);
+[ichildgrd_p,jchildgrd_p]=meshgrid(ipchild,jpchild);
+[ichildgrd_r,jchildgrd_r]=meshgrid(irchild,jrchild);
+[ichildgrd_u,jchildgrd_u]=meshgrid(ipchild,jrchild);
+[ichildgrd_v,jchildgrd_v]=meshgrid(irchild,jpchild);
+%
+% interpolations
+% 
+disp(' ')
+disp(' Do the interpolations...')                 
+np=netcdf(parent_ndepo);
+nc=netcdf(child_ndepo,'write');
+disp('ndepo...')
+for tindex=1:length(ndepot)
+  interpvar3d(np,nc,igrd_r,jgrd_r,ichildgrd_r,jchildgrd_r,'ndepo',mask,tindex)
+  interpvar3d(np,nc,igrd_r,jgrd_r,ichildgrd_r,jchildgrd_r,'noyndepo',mask,tindex)
+  interpvar3d(np,nc,igrd_r,jgrd_r,ichildgrd_r,jchildgrd_r,'noyndepo',mask,tindex)
+end
+result=close(np);
+result=close(nc);
+disp(' ')
+disp(' Done ')
+%
+% Make a plot
+%
+disp(' ')
+disp(' Make a plot...')
+figure(1)
+plot_nestndepo(child_ndepo,'ndepo',[1 6],1)
diff --git a/Nesting_tools/plot_nestndepo.m b/Nesting_tools/plot_nestndepo.m
new file mode 100644
index 00000000..111a655a
--- /dev/null
+++ b/Nesting_tools/plot_nestndepo.m
@@ -0,0 +1,83 @@
+function plot_nestndepo(child_ndepo,thefield,thetime,skip)
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%
+% Test the embedded ndepo forcing file.
+%
+%  Further Information:  
+%  http://www.croco-ocean.org
+%  
+%  This file is part of CROCOTOOLS
+%
+%  CROCOTOOLS is free software; you can redistribute it and/or modify
+%  it under the terms of the GNU General Public License as published
+%  by the Free Software Foundation; either version 2 of the License,
+%  or (at your option) any later version.
+%
+%  CROCOTOOLS is distributed in the hope that it will be useful, but
+%  WITHOUT ANY WARRANTY; without even the implied warranty of
+%  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+%  GNU General Public License for more details.
+%
+%  You should have received a copy of the GNU General Public License
+%  along with this program; if not, write to the Free Software
+%  Foundation, Inc., 59 Temple Place, Suite 330, Boston,
+%  MA  02111-1307  USA
+%
+%  Copyright (c) 2004-2006 by Pierrick Penven 
+%  e-mail:Pierrick.Penven@ird.fr  
+%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+npts=[0 0 0 0];
+i=0;
+for time=thetime
+  i=i+1;
+  
+  subplot(2,length(thetime),i)
+
+
+  nc=netcdf(child_ndepo);
+  parent_ndepo=nc.parent_file(:);
+  child_grd=nc.grd_file(:);
+  fieldc=nc{thefield}(time,:,:);
+  fieldname=nc{thefield}.long_name(:);
+  result=close(nc);
+
+  nc=netcdf(child_grd);
+  parent_grd=nc.parent_grid(:);
+  refinecoeff=nc{'refine_coef'}(:);
+  lonc=nc{'lon_rho'}(:);
+  latc=nc{'lat_rho'}(:);
+  mask=nc{'mask_rho'}(:);
+  result=close(nc);
+  mask(mask==0)=NaN;
+  pcolor(lonc,latc,mask.*fieldc)
+  shading flat
+  axis image
+  caxis([min(min(fieldc)) max(max(fieldc))])
+  colorbar
+  axis([min(min(lonc)) max(max(lonc)) min(min(latc)) max(max(latc))])
+  title(['\bf ',fieldname,' CHILD'])
+
+  subplot(2,length(thetime),i+length(thetime))
+  nc=netcdf(parent_ndepo);
+  field=nc{thefield}(time,:,:);
+  fieldname=nc{thefield}.long_name(:);
+  result=close(nc);
+  
+  nc=netcdf(parent_grd);
+  lon=nc{'lon_rho'}(:);
+  lat=nc{'lat_rho'}(:);
+  mask=nc{'mask_rho'}(:);
+  result=close(nc);
+  mask(mask==0)=NaN;
+  pcolor(lon,lat,mask.*field)
+  shading flat
+  axis image
+  caxis([min(min(fieldc)) max(max(fieldc))])
+  colorbar
+  axis([min(min(lonc)) max(max(lonc)) min(min(latc)) max(max(latc))])
+  title(['\bf ',fieldname,' PARENT'])
+end
+
+
+return
-- 
GitLab