Index: /issm/trunk-jpl/src/m/classes/SMBarma.py
===================================================================
--- /issm/trunk-jpl/src/m/classes/SMBarma.py	(revision 27483)
+++ /issm/trunk-jpl/src/m/classes/SMBarma.py	(revision 27484)
@@ -8,5 +8,5 @@
 
 class SMBarma(object):
-    """SMBARMA class definition
+    """SMBarma class definition
 
     Usage:
Index: /issm/trunk-jpl/src/m/classes/stochasticforcing.m
===================================================================
--- /issm/trunk-jpl/src/m/classes/stochasticforcing.m	(revision 27483)
+++ /issm/trunk-jpl/src/m/classes/stochasticforcing.m	(revision 27484)
@@ -111,14 +111,14 @@
 					end
 					if(strcmp(class(md.friction),'friction') || strcmp(class(md.friction),'frictionschoof') || strcmp(class(md.friction),'frictioncoulomb'))
-                  if(md.friction.coupling~=0 && md.friction.coupling~=1 && md.friction.coupling~=2)
-                     error('stochasticforcing field %s is only implemented for cases md.friction.coupling 0 or 1 or 2', char(field));
-                  end
+						if(md.friction.coupling~=0 && md.friction.coupling~=1 && md.friction.coupling~=2)
+							error('stochasticforcing field %s is only implemented for cases md.friction.coupling 0 or 1 or 2', char(field));
+						end
 					end
 					if(strcmp(class(md.friction),'friction'))
-                  if(any(md.friction.q==0))
-                     error('stochasticforcing field %s requires non-zero q exponent',char(field));
-                  end
-					end
-            end
+						if(any(md.friction.q==0))
+							error('stochasticforcing field %s requires non-zero q exponent',char(field));
+						end
+					end
+				end
 				%Checking for specific dimensions
 				if ~(strcmp(field,'SMBarma') || strcmp(field,'FrontalForcingsRignotarma') || strcmp(field,'BasalforcingsDeepwaterMeltingRatearma') || strcmp(field,'FrontalForcingsSubglacialDischargearma') || ((strcmp(field,'FrictionWaterPressure') && ispwHydroarma)))
@@ -173,81 +173,80 @@
 
 			%%% Check consistency between ARMA models %%%
-			if(indSMBarma~=-1 && indTFarma~=-1) %both ARMA models are used: check ARMA time step consistency
-				if(md.smb.arma_timestep~=md.frontalforcings.arma_timestep)
-					crossentries = reshape(self.covariance(1+sum(dimensions(1:indSMBarma-1)):sum(dimensions(1:indSMBarma)),1+sum(dimensions(1:indTFarma-1)):sum(dimensions(1:indTFarma))),1,[]);
-					if any(crossentries~=0)
-						error('SMBarma and FrontalForcingsRignotarma have different arma_timestep and non-zero covariance');
-					end
-				end
-			end
-			if(indSMBarma~=-1 && indBDWarma~=-1) %both ARMA models are used: check ARMA time step consistency
-				if(md.smb.arma_timestep~=md.basalforcings.arma_timestep)
-					crossentries = reshape(self.covariance(1+sum(dimensions(1:indSMBarma-1)):sum(dimensions(1:indSMBarma)),1+sum(dimensions(1:indBDWarma-1)):sum(dimensions(1:indBDWarma))),1,[]);
-					if any(crossentries~=0)
-						error('SMBarma and BasalforcingsDeepwaterMeltingRatearma have different arma_timestep and non-zero covariance');
-					end
-				end
-			end
-			if(indTFarma~=-1 && indBDWarma~=-1) %both ARMA models are used: check ARMA time step consistency
-				if(md.frontalforcings.arma_timestep~=md.basalforcings.arma_timestep)
-					crossentries = reshape(self.covariance(1+sum(dimensions(1:indTFarma-1)):sum(dimensions(1:indTFarma)),1+sum(dimensions(1:indBDWarma-1)):sum(dimensions(1:indBDWarma))),1,[]);
-					if any(crossentries~=0)
-						error('FrontalForcingsRignotarma and BasalforcingsDeepwaterMeltingRatearma have different arma_timestep and non-zero covariance');
-					end
-				end
-			end
-			if(indSMBarma~=-1 && indSdarma~=-1) %both ARMA models are used: check ARMA time step consistency
-				if(md.smb.arma_timestep~=md.frontalforcings.sd_arma_timestep)
-					crossentries = reshape(self.covariance(1+sum(dimensions(1:indSMBarma-1)):sum(dimensions(1:indSMBarma)),1+sum(dimensions(1:indSdarma-1)):sum(dimensions(1:indSdarma))),1,[]);
-					if any(crossentries~=0)
-						error('SMBarma and FrontalForcingsSubglacialDischargearma have different arma_timestep and non-zero covariance');
-					end
-				end
-			end
-			if(indSMBarma~=-1 && indPwarma~=-1) %both ARMA models are used: check ARMA time step consistency
-				if(md.smb.arma_timestep~=md.hydrology.arma_timestep)
-					crossentries = reshape(self.covariance(1+sum(dimensions(1:indSMBarma-1)):sum(dimensions(1:indSMBarma)),1+sum(dimensions(1:indPwarma-1)):sum(dimensions(1:indPwarma))),1,[]);
-					if any(crossentries~=0)
-						error('SMBarma and hydrologyarmapw have different arma_timestep and non-zero covariance');
-					end
-				end
-			end
-			if(indSdarma~=-1 && indBDWarma~=-1) %both ARMA models are used: check ARMA time step consistency
-				if(md.frontalforcings.sd_arma_timestep~=md.basalforcings.arma_timestep)
-					crossentries = reshape(self.covariance(1+sum(dimensions(1:indSdarma-1)):sum(dimensions(1:indSdarma)),1+sum(dimensions(1:indBDWarma-1)):sum(dimensions(1:indBDWarma))),1,[]);
-					if any(crossentries~=0)
-						error('FrontalForcingsSubglacialDischargearma and BasalforcingsDeepwaterMeltingRatearma have different arma_timestep and non-zero covariance');
-					end
-				end
-			end
-			if(indSdarma~=-1 && indTFarma~=-1) %both ARMA models are used: check ARMA time step consistency
-				if(md.frontalforcings.sd_arma_timestep~=md.frontalforcings.arma_timestep)
-					crossentries = reshape(self.covariance(1+sum(dimensions(1:indSdarma-1)):sum(dimensions(1:indSdarma)),1+sum(dimensions(1:indTFarma-1)):sum(dimensions(1:indTFarma))),1,[]);
-					if any(crossentries~=0)
-						error('FrontalForcingsRignotarma and FrontalForcingsSubglacialDischargearma have different arma_timestep and non-zero covariance');
-					end
-				end
-			end
-			if(indTFarma~=-1 && indPwarma~=-1) %both ARMA models are used: check ARMA time step consistency
-				if(md.frontalforcings.arma_timestep~=md.hydrology.arma_timestep)
-					crossentries = reshape(self.covariance(1+sum(dimensions(1:indTFarma-1)):sum(dimensions(1:indTFarma)),1+sum(dimensions(1:indPwarma-1)):sum(dimensions(1:indPwarma))),1,[]);
-					if any(crossentries~=0)
-						error('FrontalForcingsRignotarma and hydrologyarmapw have different arma_timestep and non-zero covariance');
-					end
-				end
-			end
-			if(indSdarma~=-1 && indPwarma~=-1) %both ARMA models are used: check ARMA time step consistency
-				if(md.frontalforcings.sd_arma_timestep~=md.hydrology.arma_timestep)
-					crossentries = reshape(self.covariance(1+sum(dimensions(1:indSdarma-1)):sum(dimensions(1:indSdarma)),1+sum(dimensions(1:indPwarma-1)):sum(dimensions(1:indPwarma))),1,[]);
-					if any(crossentries~=0)
-						error('FrontalForcingsSubglacialDischargearma and hydrologyarmapw have different arma_timestep and non-zero covariance');
-					end
-				end
-			end
-			if(indBDWarma~=-1 && indPwarma~=-1) %both ARMA models are used: check ARMA time step consistency
-				if(md.basalforcings.arma_timestep~=md.hydrology.arma_timestep)
-					crossentries = reshape(self.covariance(1+sum(dimensions(1:indBDWarma-1)):sum(dimensions(1:indBDWarma)),1+sum(dimensions(1:indPwarma-1)):sum(dimensions(1:indPwarma))),1,[]);
-					if any(crossentries~=0)
-						error('BasalforcingsDeepwaterMeltingRatearma and hydrologyarmapw have different arma_timestep and non-zero covariance');
+			if(indBDWarma~=-1)
+				if(indPwarma~=-1)
+					if(md.basalforcings.arma_timestep~=md.hydrology.arma_timestep)
+						crossentries = reshape(self.covariance(1+sum(dimensions(1:indBDWarma-1)):sum(dimensions(1:indBDWarma)),1+sum(dimensions(1:indPwarma-1)):sum(dimensions(1:indPwarma))),1,[]);
+						if any(crossentries~=0)
+							error('BasalforcingsDeepwaterMeltingRatearma and hydrologyarmapw have different arma_timestep and non-zero covariance');
+						end
+					end
+				elseif(indSdarma~=-1)
+					if(md.frontalforcings.sd_arma_timestep~=md.basalforcings.arma_timestep)
+						crossentries = reshape(self.covariance(1+sum(dimensions(1:indSdarma-1)):sum(dimensions(1:indSdarma)),1+sum(dimensions(1:indBDWarma-1)):sum(dimensions(1:indBDWarma))),1,[]);
+						if any(crossentries~=0)
+							error('FrontalForcingsSubglacialDischargearma and BasalforcingsDeepwaterMeltingRatearma have different arma_timestep and non-zero covariance');
+						end
+					end
+				elseif(indSMBarma~=-1)
+					if(md.smb.arma_timestep~=md.basalforcings.arma_timestep)
+						crossentries = reshape(self.covariance(1+sum(dimensions(1:indSMBarma-1)):sum(dimensions(1:indSMBarma)),1+sum(dimensions(1:indBDWarma-1)):sum(dimensions(1:indBDWarma))),1,[]);
+						if any(crossentries~=0)
+							error('SMBarma and BasalforcingsDeepwaterMeltingRatearma have different arma_timestep and non-zero covariance');
+						end
+					end
+				elseif(indTFarma~=-1)
+					if(md.frontalforcings.arma_timestep~=md.basalforcings.arma_timestep)
+						crossentries = reshape(self.covariance(1+sum(dimensions(1:indTFarma-1)):sum(dimensions(1:indTFarma)),1+sum(dimensions(1:indBDWarma-1)):sum(dimensions(1:indBDWarma))),1,[]);
+						if any(crossentries~=0)
+							error('FrontalForcingsRignotarma and BasalforcingsDeepwaterMeltingRatearma have different arma_timestep and non-zero covariance');
+						end
+					end
+				end
+			elseif(indPwarma~=-1)
+				if(indSdarma~=-1)
+					if(md.frontalforcings.sd_arma_timestep~=md.hydrology.arma_timestep)
+						crossentries = reshape(self.covariance(1+sum(dimensions(1:indSdarma-1)):sum(dimensions(1:indSdarma)),1+sum(dimensions(1:indPwarma-1)):sum(dimensions(1:indPwarma))),1,[]);
+						if any(crossentries~=0)
+							error('FrontalForcingsSubglacialDischargearma and hydrologyarmapw have different arma_timestep and non-zero covariance');
+						end
+					end
+				elseif(indSMBarma~=-1)
+					if(md.smb.arma_timestep~=md.hydrology.arma_timestep)
+						crossentries = reshape(self.covariance(1+sum(dimensions(1:indSMBarma-1)):sum(dimensions(1:indSMBarma)),1+sum(dimensions(1:indPwarma-1)):sum(dimensions(1:indPwarma))),1,[]);
+						if any(crossentries~=0)
+							error('SMBarma and hydrologyarmapw have different arma_timestep and non-zero covariance');
+						end
+					end
+				elseif(indTFarma~=-1)
+					if(md.frontalforcings.arma_timestep~=md.hydrology.arma_timestep)
+						crossentries = reshape(self.covariance(1+sum(dimensions(1:indTFarma-1)):sum(dimensions(1:indTFarma)),1+sum(dimensions(1:indPwarma-1)):sum(dimensions(1:indPwarma))),1,[]);
+						if any(crossentries~=0)
+							error('FrontalForcingsRignotarma and hydrologyarmapw have different arma_timestep and non-zero covariance');
+						end
+					end
+				end
+			elseif(indSdarma~=-1)
+				if(indSMBarma~=-1)
+					if(md.smb.arma_timestep~=md.frontalforcings.sd_arma_timestep)
+						crossentries = reshape(self.covariance(1+sum(dimensions(1:indSMBarma-1)):sum(dimensions(1:indSMBarma)),1+sum(dimensions(1:indSdarma-1)):sum(dimensions(1:indSdarma))),1,[]);
+						if any(crossentries~=0)
+							error('SMBarma and FrontalForcingsSubglacialDischargearma have different arma_timestep and non-zero covariance');
+						end
+					end
+				elseif(indTFarma~=-1)
+					if(md.frontalforcings.sd_arma_timestep~=md.frontalforcings.arma_timestep)
+						crossentries = reshape(self.covariance(1+sum(dimensions(1:indSdarma-1)):sum(dimensions(1:indSdarma)),1+sum(dimensions(1:indTFarma-1)):sum(dimensions(1:indTFarma))),1,[]);
+						if any(crossentries~=0)
+							error('FrontalForcingsRignotarma and FrontalForcingsSubglacialDischargearma have different arma_timestep and non-zero covariance');
+						end
+					end
+				end
+			elseif(indSMBarma~=-1)
+				if(indTFarma~=-1)
+					if(md.smb.arma_timestep~=md.frontalforcings.arma_timestep)
+						crossentries = reshape(self.covariance(1+sum(dimensions(1:indSMBarma-1)):sum(dimensions(1:indSMBarma)),1+sum(dimensions(1:indTFarma-1)):sum(dimensions(1:indTFarma))),1,[]);
+						if any(crossentries~=0)
+							error('SMBarma and FrontalForcingsRignotarma have different arma_timestep and non-zero covariance');
+						end
 					end
 				end
@@ -321,6 +320,6 @@
 					end
 					if(strcmp(field,'FrictionWaterPressure') && ispwHydroarma)
-                  dimensions(ind) = md.hydrology.num_basins;
-               end
+						dimensions(ind) = md.hydrology.num_basins;
+				   end
 					ind = ind+1;
 				end
Index: /issm/trunk-jpl/src/m/classes/stochasticforcing.py
===================================================================
--- /issm/trunk-jpl/src/m/classes/stochasticforcing.py	(revision 27483)
+++ /issm/trunk-jpl/src/m/classes/stochasticforcing.py	(revision 27484)
@@ -7,5 +7,5 @@
 
 class stochasticforcing(object):
-    """STOCHASTICFORCING class definition
+    """stochasticforcing class definition
 
     Usage:
@@ -65,5 +65,5 @@
 
         num_fields = len(self.fields)
-        if(self.stochastictimestep==0):
+        if self.stochastictimestep == 0:
             md.stochasticforcing.stochastictimestep = md.timestepping.time_step #by default: stochastictimestep set to ISSM time step
             print('      stochasticforcing.stocahstictimestep not specified: set to md.timestepping.time_step')
@@ -79,5 +79,5 @@
         structstoch = self.structstochforcing()
         # Check if hydrologyarmapw is used
-        if((type(md.hydrology).__name__ == 'hydrologyarmapw') and md.transient.hydrology==1):
+        if type(md.hydrology).__name__ == 'hydrologyarmapw' and md.transient.ishydrology:
             ispwHydroarma = 1
         else:
@@ -87,29 +87,29 @@
             if 'SMBarma' in field:
                 mdname = structstoch[field]
-                if (type(md.smb).__name__ != mdname):
+                if type(md.smb).__name__ != mdname:
                     raise TypeError('md.smb does not agree with stochasticforcing field {}'.format(field))
             if 'SMBforcing' in field:
                 mdname = structstoch[field]
-                if (type(md.smb).__name__ != mdname):
+                if type(md.smb).__name__ != mdname:
                     raise TypeError('md.smb does not agree with stochasticforcing field {}'.format(field))
             if 'FrontalForcings' in field:
                 mdname = structstoch[field]
-                if (type(md.frontalforcings).__name__ != mdname):
+                if type(md.frontalforcings).__name__ != mdname:
                     raise TypeError('md.frontalforcings does not agree with stochasticforcing field {}'.format(field))
             if 'Calving' in field:
                 mdname = structstoch[field]
-                if (type(md.calving).__name__ != mdname):
+                if type(md.calving).__name__ != mdname:
                     raise TypeError('md.calving does not agree with stochasticforcing field {}'.format(field))
             if 'BasalforcingsFloatingice' in field:
                 mdname = structstoch[field]
-                if (type(md.basalforcings).__name__ != mdname):
+                if type(md.basalforcings).__name__ != mdname:
                     raise TypeError('md.basalforcings does not agree with stochasticforcing field {}'.format(field))
             if 'BasalforcingsSpatialDeepwaterMeltingRate' in field:
                 mdname = structstoch[field]
-                if (type(md.basalforcings).__name__ != mdname):
+                if type(md.basalforcings).__name__ != mdname:
                     raise TypeError('md.basalforcings does not agree with stochasticforcing field {}'.format(field))
             if 'BasalforcingsDeepwaterMeltingRatearma' in field:
                 mdname = structstoch[field]
-                if (type(md.basalforcings).__name__ != mdname):
+                if type(md.basalforcings).__name__ != mdname:
                     raise TypeError('md.basalforcings does not agree with stochasticforcing field {}'.format(field))
             if 'WaterPressure' in field:
@@ -118,23 +118,23 @@
                 found   = 0
                 for ii in range(len(mdnames)):
-                    if (type(md.friction).__name__ == mdnames[ii]):
+                    if type(md.friction).__name__ == mdnames[ii]:
                         found = 1
-                if(found==0):
+                if not found:
                     raise TypeError('md.friction does not agree with stochasticforcing field {}'.format(field))
                 #if (type(md.friction).__name__ != mdname):
                 #    raise TypeError('md.friction does not agree with stochasticforcing field {}'.format(field))
-                if (type(md.friction).__name__=='friction' or type(md.friction).__name__=='frictionschoof' or type(md.friction).__name__=='frictioncoulomb'):
+                if type(md.friction).__name__ == 'friction' or type(md.friction).__name__ == 'frictionschoof' or type(md.friction).__name__=='frictioncoulomb':
                     if md.friction.coupling not in[0, 1, 2]:
                         raise TypeError('stochasticforcing field {} is only implemented for cases md.friction.coupling 0 or 1 or 2'.format(field))
-                if (type(md.friction).__name__=='friction'):
+                if type(md.friction).__name__ == 'friction':
                     if (np.any(md.friction.q == 0)):
                         raise TypeError('stochasticforcing field {} requires non-zero q exponent'.format(field))
 
             # Checking for specific dimensions
-            if ((field not in['SMBarma', 'FrontalForcingsRignotarma','BasalforcingsDeepwaterMeltingRatearma']) and ((field=='FrictionWaterPressure' and ispwHydroarma==True)==False)):
+            if field not in ['SMBarma', 'FrontalForcingsRignotarma','BasalforcingsDeepwaterMeltingRatearma']  and not (field == 'FrictionWaterPressure' and ispwHydroarma == True):
                 checkdefaults = True  # field with non-specific dimensionality
 
         # Retrieve sum of all the field dimensionalities
-        dimensions = self.defaultdimension*np.ones((num_fields))
+        dimensions = self.defaultdimension * np.ones((num_fields))
         indSMBarma   = -1  # About to check for index of SMBarma
         indTFarma    = -1  # About to check for index of FrontalForcingsRignotarma
@@ -142,72 +142,75 @@
         indBDWarma   = -1  # About to check for index of BasalforcingsDeepwaterMeltingRatearma
         indPwarma    = -1  # About to check for index of hydrologyarmapw
-        if ('SMBarma' in self.fields):
+        if 'SMBarma' in self.fields:
             indSMBarma = self.fields.index('SMBarma')  # Index of SMBarma, now check for consistency with other timesteps
             dimensions[indSMBarma] = md.smb.num_basins
             if(md.smb.arma_timestep<self.stochastictimestep):
                 raise TypeError('SMBarma cannot have a timestep shorter than stochastictimestep')
-        if ('FrontalForcingsRignotarma' in self.fields):
+        if 'FrontalForcingsRignotarma' in self.fields:
             indTFarma = self.fields.index('FrontalForcingsRignotarma')  # Index of TFarma, now check for consistency with other timesteps
             dimensions[indTFarma] = md.frontalforcings.num_basins
-            if(md.frontalforcings.arma_timestep<self.stochastictimestep):
+            if md.frontalforcings.arma_timestep < self.stochastictimestep:
                 raise TypeError('FrontalForcingsRignotarma cannot have a timestep shorter than stochastictimestep')
-        if ('FrontalForcingsSubglacialDischargearma' in self.fields):
+        if 'FrontalForcingsSubglacialDischargearma' in self.fields:
             indSdarma = self.fields.index('FrontalForcingsSubglacialDischargearma')  # Index of Sdarma, now check for consistency with other timesteps
             dimensions[indSdarma] = md.frontalforcings.num_basins
-            if(md.frontalforcings.sd_arma_timestep<self.stochastictimestep):
+            if md.frontalforcings.sd_arma_timestep < self.stochastictimestep:
                 raise TypeError('FrontalForcingsSubglacialDischargearma cannot have a timestep shorter than stochastictimestep')
-        if ('BasalforcingsDeepwaterMeltingRatearma' in self.fields):
+        if 'BasalforcingsDeepwaterMeltingRatearma' in self.fields:
             indBDWarma = self.fields.index('BasalforcingsDeepwaterMeltingRatearma')  # Index of BDWarma, now check for consistency with other timesteps
             dimensions[indTFarma] = md.basalforcings.num_basins
-            if(md.basalforcings.arma_timestep<self.stochastictimestep):
+            if md.basalforcings.arma_timestep < self.stochastictimestep:
                 raise TypeError('BasalforcingsDeepwaterMeltingRatearma cannot have a timestep shorter than stochastictimestep')
-        if (('FrictionWaterPressure' in self.fields) and ispwHydroarma):
+        if 'FrictionWaterPressure' in self.fields and ispwHydroarma:
             indPwarma = self.fields.index('FrictionWaterPressure')  # Index of Pwarma, now check for consistency with other timesteps
             dimensions[indPwarma] = md.hydrology.num_basins
-            if(md.hydrology.arma_timestep<self.stochastictimestep):
+            if md.hydrology.arma_timestep < self.stochastictimestep:
                 raise TypeError('hydrologyarmapw cannot have a timestep shorter than stochastictimestep')
         size_tot = np.sum(dimensions)
 
-        if (indSMBarma != -1 and indTFarma != -1):  # Both ARMA models are used: check ARMA time step consistency
-            covsum = self.covariance[np.sum(dimensions[0:indSMBarma]).astype(int):np.sum(dimensions[0:indSMBarma + 1]).astype(int), np.sum(dimensions[0:indTFarma]).astype(int):np.sum(dimensions[0:indTFarma + 1]).astype(int)]
-            if((md.smb.arma_timestep != md.frontalforcings.arma_timestep) and np.any(covsum != 0)):
-                raise IOError('SMBarma and FrontalForcingsRignotarma have different arma_timestep and non-zero covariance')
-        if (indSMBarma != -1 and indBDWarma != -1):  # Both ARMA models are used: check ARMA time step consistency
-            covsum = self.covariance[np.sum(dimensions[0:indSMBarma]).astype(int):np.sum(dimensions[0:indSMBarma + 1]).astype(int), np.sum(dimensions[0:indBDWarma]).astype(int):np.sum(dimensions[0:indBDWarma + 1]).astype(int)]
-            if((md.smb.arma_timestep != md.basalforcings.arma_timestep) and np.any(covsum != 0)):
-                raise IOError('SMBarma and BasalforcingsDeepwaterMeltingRatearma have different arma_timestep and non-zero covariance')
-        if (indTFarma != -1 and indBDWarma != -1):  # Both ARMA models are used: check ARMA time step consistency
-            covsum = self.covariance[np.sum(dimensions[0:indTFarma]).astype(int):np.sum(dimensions[0:indTFarma + 1]).astype(int), np.sum(dimensions[0:indBDWarma]).astype(int):np.sum(dimensions[0:indBDWarma + 1]).astype(int)]
-            if((md.frontalforcings.arma_timestep != md.basalforcings.arma_timestep) and np.any(covsum != 0)):
-                raise IOError('FrontalForcingsRignotarma and BasalforcingsDeepwaterMeltingRatearma have different arma_timestep and non-zero covariance')
-        if (indSMBarma != -1 and indSdarma != -1):  # Both ARMA models are used: check ARMA time step consistency
-            covsum = self.covariance[np.sum(dimensions[0:indSMBarma]).astype(int):np.sum(dimensions[0:indSMBarma + 1]).astype(int), np.sum(dimensions[0:indSdarma]).astype(int):np.sum(dimensions[0:indSdarma + 1]).astype(int)]
-            if((md.smb.arma_timestep != md.frontalforcings.sd_arma_timestep) and np.any(covsum != 0)):
-                raise IOError('SMBarma and FrontalForcingsSubglacialDischargearma have different arma_timestep and non-zero covariance')
-        if (indSdarma != -1 and indBDWarma != -1):  # Both ARMA models are used: check ARMA time step consistency
-            covsum = self.covariance[np.sum(dimensions[0:indSdarma]).astype(int):np.sum(dimensions[0:indSdarma + 1]).astype(int), np.sum(dimensions[0:indBDWarma]).astype(int):np.sum(dimensions[0:indBDWarma + 1]).astype(int)]
-            if((md.frontalforcings.sd_arma_timestep != md.basalforcings.arma_timestep) and np.any(covsum != 0)):
-                raise IOError('FrontalForcingsSubglacialDischargearma and BasalforcingsDeepwaterMeltingRatearma have different arma_timestep and non-zero covariance')
-        if (indSdarma != -1 and indTFarma != -1):  # Both ARMA models are used: check ARMA time step consistency
-            covsum = self.covariance[np.sum(dimensions[0:indSdarma]).astype(int):np.sum(dimensions[0:indSdarma + 1]).astype(int), np.sum(dimensions[0:indTFarma]).astype(int):np.sum(dimensions[0:indTFarma + 1]).astype(int)]
-            if((md.frontalforcings.sd_arma_timestep != md.frontalforcings.arma_timestep) and np.any(covsum != 0)):
-                raise IOError('FrontalForcingsSubglacialDischargearma and FrontalForcingsRignotarma have different arma_timestep and non-zero covariance')
-        if (indSMBarma != -1 and indPwarma != -1):  # Both ARMA models are used: check ARMA time step consistency
-            covsum = self.covariance[np.sum(dimensions[0:indSMBarma]).astype(int):np.sum(dimensions[0:indSMBarma + 1]).astype(int), np.sum(dimensions[0:indPwarma]).astype(int):np.sum(dimensions[0:indPwarma + 1]).astype(int)]
-            if((md.smb.arma_timestep != md.hydrology.arma_timestep) and np.any(covsum != 0)):
-                raise IOError('SMBarma and hydrologyarmapw have different arma_timestep and non-zero covariance')
-        if (indTFarma != -1 and indPwarma != -1):  # Both ARMA models are used: check ARMA time step consistency
-            covsum = self.covariance[np.sum(dimensions[0:indTFarma]).astype(int):np.sum(dimensions[0:indTFarma + 1]).astype(int), np.sum(dimensions[0:indPwarma]).astype(int):np.sum(dimensions[0:indPwarma + 1]).astype(int)]
-            if((md.frontalforcings.arma_timestep != md.hydrology.arma_timestep) and np.any(covsum != 0)):
-                raise IOError('FrontalForcingsRignotarma and hydrologyarmapw have different arma_timestep and non-zero covariance')
-        if (indSdarma != -1 and indPwarma != -1):  # Both ARMA models are used: check ARMA time step consistency
-            covsum = self.covariance[np.sum(dimensions[0:indSdarma]).astype(int):np.sum(dimensions[0:indSdarma + 1]).astype(int), np.sum(dimensions[0:indPwarma]).astype(int):np.sum(dimensions[0:indPwarma + 1]).astype(int)]
-            if((md.frontalforcings.sd_arma_timestep != md.hydrology.arma_timestep) and np.any(covsum != 0)):
-                raise IOError('FrontalForingsSubglacialDischargearma and hydrologyarmapw have different arma_timestep and non-zero covariance')
-        if (indBDWarma != -1 and indPwarma != -1):  # Both ARMA models are used: check ARMA time step consistency
-            covsum = self.covariance[np.sum(dimensions[0:indBDWarma]).astype(int):np.sum(dimensions[0:indBDWarma + 1]).astype(int), np.sum(dimensions[0:indPwarma]).astype(int):np.sum(dimensions[0:indPwarma + 1]).astype(int)]
-            if((md.smb.arma_timestep != md.hydrology.arma_timestep) and np.any(covsum != 0)):
-                raise IOError('BasalforcingsDeepwaterMeltingRatarma and hydrologyarmapw have different arma_timestep and non-zero covariance')
- 
+        if indBDWarma != -1:
+            if indPwarma != -1:
+                covsum = self.covariance[np.sum(dimensions[0:indBDWarma]).astype(int):np.sum(dimensions[0:indBDWarma + 1]).astype(int), np.sum(dimensions[0:indPwarma]).astype(int):np.sum(dimensions[0:indPwarma + 1]).astype(int)]
+                if md.smb.arma_timestep != md.hydrology.arma_timestep and np.any(covsum != 0):
+                    raise IOError('BasalforcingsDeepwaterMeltingRatarma and hydrologyarmapw have different arma_timestep and non-zero covariance')
+            elif indSdarma != -1:
+                covsum = self.covariance[np.sum(dimensions[0:indSdarma]).astype(int):np.sum(dimensions[0:indSdarma + 1]).astype(int), np.sum(dimensions[0:indBDWarma]).astype(int):np.sum(dimensions[0:indBDWarma + 1]).astype(int)]
+                if md.frontalforcings.sd_arma_timestep != md.basalforcings.arma_timestep and np.any(covsum != 0):
+                    raise IOError('FrontalForcingsSubglacialDischargearma and BasalforcingsDeepwaterMeltingRatearma have different arma_timestep and non-zero covariance')
+            elif indSMBarma != -1:
+                covsum = self.covariance[np.sum(dimensions[0:indSMBarma]).astype(int):np.sum(dimensions[0:indSMBarma + 1]).astype(int), np.sum(dimensions[0:indBDWarma]).astype(int):np.sum(dimensions[0:indBDWarma + 1]).astype(int)]
+                if md.smb.arma_timestep != md.basalforcings.arma_timestep and np.any(covsum != 0):
+                    raise IOError('SMBarma and BasalforcingsDeepwaterMeltingRatearma have different arma_timestep and non-zero covariance')
+            elif indTFarma != -1:
+                covsum = self.covariance[np.sum(dimensions[0:indTFarma]).astype(int):np.sum(dimensions[0:indTFarma + 1]).astype(int), np.sum(dimensions[0:indBDWarma]).astype(int):np.sum(dimensions[0:indBDWarma + 1]).astype(int)]
+                if md.frontalforcings.arma_timestep != md.basalforcings.arma_timestep and np.any(covsum != 0):
+                    raise IOError('FrontalForcingsRignotarma and BasalforcingsDeepwaterMeltingRatearma have different arma_timestep and non-zero covariance')
+        elif indPwarma != -1:
+            if indSdarma != -1:
+                covsum = self.covariance[np.sum(dimensions[0:indSdarma]).astype(int):np.sum(dimensions[0:indSdarma + 1]).astype(int), np.sum(dimensions[0:indPwarma]).astype(int):np.sum(dimensions[0:indPwarma + 1]).astype(int)]
+                if md.frontalforcings.sd_arma_timestep != md.hydrology.arma_timestep and np.any(covsum != 0):
+                    raise IOError('FrontalForingsSubglacialDischargearma and hydrologyarmapw have different arma_timestep and non-zero covariance')
+            elif indSMBarma != -1:
+                covsum = self.covariance[np.sum(dimensions[0:indSMBarma]).astype(int):np.sum(dimensions[0:indSMBarma + 1]).astype(int), np.sum(dimensions[0:indPwarma]).astype(int):np.sum(dimensions[0:indPwarma + 1]).astype(int)]
+                if md.smb.arma_timestep != md.hydrology.arma_timestep and np.any(covsum != 0):
+                    raise IOError('SMBarma and hydrologyarmapw have different arma_timestep and non-zero covariance')
+            elif indTFarma != -1:
+                covsum = self.covariance[np.sum(dimensions[0:indTFarma]).astype(int):np.sum(dimensions[0:indTFarma + 1]).astype(int), np.sum(dimensions[0:indPwarma]).astype(int):np.sum(dimensions[0:indPwarma + 1]).astype(int)]
+                if md.frontalforcings.arma_timestep != md.hydrology.arma_timestep and np.any(covsum != 0):
+                    raise IOError('FrontalForcingsRignotarma and hydrologyarmapw have different arma_timestep and non-zero covariance')
+        elif indSdarma != -1:
+            if indSMBarma != -1:
+                covsum = self.covariance[np.sum(dimensions[0:indSMBarma]).astype(int):np.sum(dimensions[0:indSMBarma + 1]).astype(int), np.sum(dimensions[0:indSdarma]).astype(int):np.sum(dimensions[0:indSdarma + 1]).astype(int)]
+                if md.smb.arma_timestep != md.frontalforcings.sd_arma_timestep and np.any(covsum != 0):
+                    raise IOError('SMBarma and FrontalForcingsSubglacialDischargearma have different arma_timestep and non-zero covariance')
+            elif indTFarma != -1:
+                covsum = self.covariance[np.sum(dimensions[0:indSdarma]).astype(int):np.sum(dimensions[0:indSdarma + 1]).astype(int), np.sum(dimensions[0:indTFarma]).astype(int):np.sum(dimensions[0:indTFarma + 1]).astype(int)]
+                if md.frontalforcings.sd_arma_timestep != md.frontalforcings.arma_timestep and np.any(covsum != 0):
+                    raise IOError('FrontalForcingsSubglacialDischargearma and FrontalForcingsRignotarma have different arma_timestep and non-zero covariance')
+        elif indSMBarma != -1:
+            if indTFarma != -1:
+                covsum = self.covariance[np.sum(dimensions[0:indSMBarma]).astype(int):np.sum(dimensions[0:indSMBarma + 1]).astype(int), np.sum(dimensions[0:indTFarma]).astype(int):np.sum(dimensions[0:indTFarma + 1]).astype(int)]
+                if md.smb.arma_timestep != md.frontalforcings.arma_timestep and np.any(covsum != 0):
+                    raise IOError('SMBarma and FrontalForcingsRignotarma have different arma_timestep and non-zero covariance')
 
         md = checkfield(md, 'fieldname', 'stochasticforcing.isstochasticforcing', 'values', [0, 1])
@@ -216,5 +219,5 @@
         md = checkfield(md, 'fieldname', 'stochasticforcing.stochastictimestep', 'NaN', 1,'Inf', 1, '>=', md.timestepping.time_step)
         md = checkfield(md, 'fieldname', 'stochasticforcing.randomflag', 'numel', [1], 'values', [0, 1])
-        if (checkdefaults):
+        if checkdefaults:
             md = checkfield(md, 'fieldname', 'stochasticforcing.defaultdimension', 'numel', 1, 'NaN', 1, 'Inf', 1, '>', 0)
             md = checkfield(md, 'fieldname', 'stochasticforcing.default_id', 'Inf', 1, 'NaN', 1, '>=', 0, '<=', self.defaultdimension, 'size', [md.mesh.numberofelements])
@@ -235,8 +238,8 @@
         else:
             num_fields = len(self.fields)
-            if(self.stochastictimestep==0):
+            if self.stochastictimestep == 0:
                 md.stochasticforcing.stochastictimestep = md.timestepping.time_step #by default: stochastictimestep set to ISSM time step
             # Check if hydroarmapw is used
-            if((type(md.hydrology).__name__ == 'hydrologyarmapw') and md.transient.hydrology==1):
+            if type(md.hydrology).__name__ == 'hydrologyarmapw' and md.transient.ishydrology:
                 ispwHydroarma = 1
             else:
@@ -247,17 +250,17 @@
             for ind, field in enumerate(self.fields):
                 # Checking for specific dimensions
-                if (field == 'SMBarma'):
+                if field == 'SMBarma':
                     dimensions[ind] = md.smb.num_basins
-                if (field == 'FrontalForcingsRignotarma'):
+                elif field == 'FrontalForcingsRignotarma':
                     dimensions[ind] = md.frontalforcings.num_basins
-                if (field == 'FrontalForcingsSubglacialDischargearma'):
+                elif field == 'FrontalForcingsSubglacialDischargearma':
                     dimensions[ind] = md.frontalforcings.num_basins
-                if (field == 'BasalforcingsDeepwaterMeltingRatearma'):
+                elif field == 'BasalforcingsDeepwaterMeltingRatearma':
                     dimensions[ind] = md.basalforcings.num_basins
-                if (field == 'FrictionWaterPressure' and ispwHydroarma):
+                elif field == 'FrictionWaterPressure' and ispwHydroarma:
                     dimensions[ind] = md.hydrology.num_basins
 
             # Scaling covariance matrix (scale column-by-column and row-by-row)
-            scaledfields = ['BasalforcingsDeepwaterMeltingRatearma','BasalforcingsSpatialDeepwaterMeltingRate','DefaultCalving', 'FloatingMeltRate', 'SMBarma', 'SMBforcing']  # list of fields that need scaling * 1/yts
+            scaledfields = ['BasalforcingsDeepwaterMeltingRatearma','BasalforcingsSpatialDeepwaterMeltingRate','DefaultCalving', 'FloatingMeltRate', 'SMBarma', 'SMBforcing']  # list of fields that need scaling * 1 / yts
             tempcovariance = np.copy(self.covariance)
             for i in range(num_fields):