--- MITgcm_contrib/gmaze_pv/compute_MLD.m	2006/10/06 18:54:35	1.1
+++ MITgcm_contrib/gmaze_pv/compute_MLD.m	2007/09/19 15:37:38	1.6
@@ -1,5 +1,5 @@
 %
-% [] = compute_MLD(SNAPSHOT)
+% [MLD] = compute_MLD(SNAPSHOT)
 %
 % Here we compute the Mixed Layer Depth as:
 % MLD = min depth for which : ST > ST(SSS,SST-0.8,p0)  
@@ -9,7 +9,7 @@
 %  SST the Sea Surface Temperature (oC)
 %  SSS the Sea Surface Salinity (PSU-35)
 %  p0  the Sea Level Pressure (mb)
-%  EKL is the Ekman layer depth (m)
+%  EKL is the Ekman layer depth (m, positive)
 %
 % Files names are:
 % INPUT:
@@ -28,7 +28,7 @@
 % 09/20/06
 % gmaze@mit.edu
 
-function compute_MLD(snapshot)
+function varargout = compute_MLD(snapshot)
 
 global sla toshow
 global netcdf_suff netcdf_domain
@@ -86,6 +86,16 @@
 
 MLD(isnan(squeeze(ST(1,:,:)))) = NaN;
 
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% Ensure we have the right sign (positive)
+mm = nanmean(nanmean(MLD,1));
+if mm <= 0
+  MLD = -MLD;
+end
+
+
+
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%
 % Record
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%
@@ -145,5 +155,14 @@
 nc{ncid}(:,:,:)        = MLD;
 
 nc=close(nc);
+close(ncST);
+close(ncS);
+close(ncT);
 
 
+% Output:
+output = struct('MLD',MLD,'lat',STlat,'lon',STlon);
+switch nargout
+ case 1
+  varargout(1) = {output};
+end