Updates for paper

README.md

@@ -23,6 +23,8 @@ You can then activate the oml environment and follow the installation instructio
 
 ### Optional additional packages to use all scripts and notebooks:  
 
+*Warning* for some reason conda does not appear to successfully build these packages after gsw is installed.  Building these packages before installing gsw does seem to work.  This error should be better understood.  
+
 ```
 conda install matplotlib jupyter netcdf4 xarray ipykernel  
 ```

oceanmixedlayers/column.py

@@ -15,9 +15,9 @@ class column():
     def __init__(self,
                  kind='idealized',
                  idealized_type='linear',T0=0,dTdz=0,S0=0,dSdz=0,Sb0=0,Tb0=0,Smin=0,Tmin=0,mixedfrac=0.5,
-                 Boussinesq=False,Compressible=False,EOS='Full',rho0=1025,grav=9.81,
+                 Boussinesq=False,Compressible=False,EOS='gsw',rho0=1025,grav=9.81,
                  nz=1,Dpt=1,
-                 ArgoPath='/net3/bgr/Datasets/Argo/202011-ArgoData/dac/aoml/',ArgoID=3900660,NP=0,MaxP=1.e8,
+                 ArgoPath='',ArgoID=3900660,NP=0,MaxP=1.e8,
                  zc=[],dZ=[],T=[],S=[],
                  Debug=False
                 ):
@@ -177,7 +177,7 @@ def GetPressure(self):
             self.dp = self.pi[:-1]-self.pi[1:]
 
     def GetRho(self):
-        if self.EOS=='Full':
+        if self.EOS=='gsw':
             self.rho=gsw.density.rho(self.S,self.T,self.pc/10000.)
             self.prho=gsw.density.rho(self.S,self.T,0.)
         if self.EOS=='Linear':

oceanmixedlayers/energy.py

@@ -48,7 +48,7 @@ def __init__(self,
                  Debug=False,
                  #Constants that can be set
                  grav=9.81,rho0=1025.,
-                 eqstate='Full'):
+                 eqstate='gsw'):
         self.grav = grav
         self.rho0 = rho0
         self.eqstate = eqstate
@@ -70,7 +70,7 @@ def __init__(self,
             Pc = Vc
             dP = dV
             #Sets Rho_i
-            if eqstate=='Full':
+            if eqstate=='gsw':
                 Rho_i = gsw.density.rho(S_layer,T_layer,Pc/1.e4)
             elif eqstate=='Linear':
                 Rho_i = self.rho0+(S_layer-35)*0.8-(T_layer-10)*0.2
@@ -96,7 +96,7 @@ def __init__(self,
             Pi = np.array([0.,]+list(np.cumsum(dZ_i*self.grav*self.rho0)))
             Pc = 0.5*(Pi[1:]+Pi[:-1])
             #Compute Rho_i assuming hydrostatic pressure
-            if eqstate=='Full':
+            if eqstate=='gsw':
                 Rho_i = gsw.density.rho(S_layer,T_layer,Pc/1.e4)
             elif eqstate=='Linear':
                 Rho_i = self.rho0+(S_layer-35)*0.8-(T_layer-10)*0.2
@@ -107,7 +107,7 @@ def __init__(self,
                     #How many iterations?
                     Pi = np.array([0.,]+list(np.cumsum(dZ_i*self.grav*Rho_i)))
                     Pc = 0.5*(Pi[1:]+Pi[:-1])
-                    if eqstate=='Full':
+                    if eqstate=='gsw':
                         Rho_0 = gsw.density.rho(S_layer,T_layer,Pc/1.e4)
                     elif eqstate=='Linear':
                         Rho_0 = self.rho0+(S_layer-35)*0.8-(T_layer-10)*0.2
@@ -207,7 +207,7 @@ def compute_MLD(self,T_i,S_i,Rho_i,Pc,dP,Zc_i,dZ_i,energy,Debug=False):
             S_x[:z+1,ACTIVE],_ = MixLayers(S_i[:,ACTIVE],
                                            -dP[:,ACTIVE],
                                            z+1)
-            if self.eqstate=='Full':
+            if self.eqstate=='gsw':
                 Rho_x = gsw.density.rho(S_x[:z+1,ACTIVE],
                                         T_x[:z+1,ACTIVE],
                                         Pc[:z+1,ACTIVE]/1.e4)
@@ -264,7 +264,7 @@ def compute_MLD(self,T_i,S_i,Rho_i,Pc,dP,Zc_i,dZ_i,energy,Debug=False):
 
 
                 #Recompute the reduced bottom layer thickness for the initial state
-                if self.eqstate=='Full':
+                if self.eqstate=='gsw':
                     Rho_i_bot = gsw.density.rho(S_i[z,FINAL],T_i[z,FINAL],PC[z,...]/1.e4)
                 elif self.eqstate=='Linear':
                     Rho_i_bot = self.rho0+(S_i[z,FINAL]-35)*0.8-(T_i[z,FINAL]-10)*0.2
@@ -278,7 +278,7 @@ def compute_MLD(self,T_i,S_i,Rho_i,Pc,dP,Zc_i,dZ_i,energy,Debug=False):
                 S_x[:z+1,FINAL],MLDp[FINAL] = MixLayers(S_i[:z+1,FINAL],
                                                  -DP,
                                                  z+1)
-                if self.eqstate=='Full':
+                if self.eqstate=='gsw':
                     Rho_x = gsw.density.rho(S_x[:z+1,FINAL],
                                             T_x[:z+1,FINAL],
                                             PC/1.e4)

oceanmixedlayers/gradient.py

@@ -291,8 +291,9 @@ def nd_linefit(coordinate,value):
     covariance[nzero] = np.nanmean((coordinate[:,nzero]-coor_mean[nzero])
                                    *(value[:,nzero]-val_mean[nzero]),axis=0)
 
-    correlation = np.zeros(np.shape(val_mean))
-    correlation[nzero] = covariance[nzero]/(val_std[nzero]*coor_std[nzero])
+    #This was computed but I don't think it is used so...
+    #correlation = np.zeros(np.shape(val_mean))
+    #correlation[nzero] = covariance[nzero]/(val_std[nzero]*coor_std[nzero])
 
     slope = np.zeros(np.shape(val_mean))
     slope[nzero] = covariance[nzero]/(coor_std[nzero]**2)