Index: /issm/trunk-jpl/src/m/contrib/larour/UQ/sampling.m
===================================================================
--- /issm/trunk-jpl/src/m/contrib/larour/UQ/sampling.m	(revision 25326)
+++ /issm/trunk-jpl/src/m/contrib/larour/UQ/sampling.m	(revision 25326)
@@ -0,0 +1,29 @@
+function s=drawfrompdf(x,y,nsamples,mode)
+% Samples a 1-D pdf histogram
+% x: variable to be sampled
+% y: pdf histogram of x, best if y(1)=0 in continuous sample mode
+% nsamples: how many samples must be generated
+% mode: 'continuous' or 'discrete'. 'discrete' will generate only values provided in x, 'continuous' will interpolate such that values provided be < x(end), or in [x(1) x(end)] if y(1)=0 
+ y=[cumsum(y/sum(y))]; %compute the CDF
+ s=zeros(nsamples,1);
+
+if y(1)~=0 & strcmpi(mode,'continuous')
+    disp('WARNING: y(1)~=0, values < x(1) might be generated extrapolating y(2)-y(1)')
+end
+
+ for i=1:nsamples
+     r=rand(1); %uniform [0 1] sampling
+    j=2;    
+    while r>y(j)
+        j=j+1;
+    end
+    if strcmpi(mode,'continuous')
+    d=(r-y(j-1))/(y(j)-y(j-1));
+    s(i)=x(j-1)*(1.d0-d)+x(j)*d; %linear interpolation of the cdf to map r into x
+    elseif strcmpi(mode, 'discrete')
+        s(i)=x(j);
+    else 
+        error(' ''mode'' must be set to ''discrete'' or ''continuous'' ')
+    end
+ end
+ 
Index: /issm/trunk-jpl/src/m/contrib/larour/equiprobable_histogram_uncertain.m
===================================================================
--- /issm/trunk-jpl/src/m/contrib/larour/equiprobable_histogram_uncertain.m	(revision 25325)
+++ /issm/trunk-jpl/src/m/contrib/larour/equiprobable_histogram_uncertain.m	(revision 25326)
@@ -1,6 +1,6 @@
-function [ids,probs,n]=equiprobable_histogram_uncertain(numberids)
+function [abscissas,counts,pair_per_variable]=equiprobable_histogram_uncertain(numberabscissas)
 
-	n=numberids+1;
-	ids=(1:(numberids+1))';
-	p=ones(numberids,1)/numberids;
-	probs=[p;0]; probs=probs/sum(probs);
+	pair_per_variable=numberabscissas+1;
+	abscissas=(1:(numberabscissas+1))';
+	p=ones(numberabscissas,1)/numberabscissas;
+	counts=[p;0]; counts=counts/sum(counts);
Index: /issm/trunk-jpl/src/m/contrib/larour/glaciermip.m
===================================================================
--- /issm/trunk-jpl/src/m/contrib/larour/glaciermip.m	(revision 25325)
+++ /issm/trunk-jpl/src/m/contrib/larour/glaciermip.m	(revision 25326)
@@ -178,4 +178,13 @@
 			end
 		end % }}}
+		function part=partition(self,md,rgi2mesh,part,value) % {{{
+
+			for i=1:size(rgi2mesh,2),
+				dh=rgi2mesh(:,i);
+				pos=find(dh);
+				part(pos)=value;
+			end
+
+		end  % }}}
 	end 
 end
Index: /issm/trunk-jpl/src/m/contrib/larour/oceanmip.m
===================================================================
--- /issm/trunk-jpl/src/m/contrib/larour/oceanmip.m	(revision 25325)
+++ /issm/trunk-jpl/src/m/contrib/larour/oceanmip.m	(revision 25326)
@@ -207,32 +207,49 @@
 		end
 		%}}}
-	function arrays= dbottompressures(self) % {{{
+	function arrays= dbottompressures(self,varargin) % {{{
+		options=pairoptions(varargin{:});
+		units=getfieldvalue(options,'units','mm/yr');
+		
 		arrays=cell(self.n,1);
 		for i=1:self.n,
-			pbo=self.mesh_pbo{i}; pbo=pbo/1000; %in meters
+			pbo=self.mesh_pbo{i}; 
 			dpbo=diff(pbo,1,2);
-			time=self.time{i};
-			time=time(1:end-1);
-			array=[dpbo;dtime];
+			if strcmpi(units,'m/yr'),
+				dpbo=dpbo/1000;
+			end
+			time=self.time{i}; time=time(1:end-1);
+			array=[dpbo;time'];
 			arrays{i}=array;
 		end
 	end %}}}
-	function arrays= dzoss(self) % {{{
+	function arrays= dzoss(self,varargin) % {{{
+		options=pairoptions(varargin{:});
+		units=getfieldvalue(options,'units','mm/yr');
+		
 		arrays=cell(self.n,1);
 		for i=1:self.n,
-			zos=self.mesh_zos{i}; zos=zos/1000; %in meters
+			zos=self.mesh_zos{i}; 
 			dzos=diff(zos,1,2);
+			if strcmpi(units,'m/yr'),
+				dzos=dzos/1000;
+			end
 			time=self.time{i};time=time(1:end-1);
-			array=[dzos;time];
+			array=[dzos;time'];
 			arrays{i}=array;
 		end
 	end	%}}}
-	function arrays= dzostogass(self) % {{{
+	function arrays= dzostogass(self,varargin) % {{{
+		options=pairoptions(varargin{:});
+		units=getfieldvalue(options,'units','mm/yr');
+
 		arrays=cell(self.n,1);
 		for i=1:self.n,
-			zostoga=self.zostoga{i}; zos=zos/1000; %in meters
-			dzostoga=diff(dzostoga);
+			zostoga=self.zostoga{i}; 
+			dzostoga=diff(zostoga);
+			if strcmpi(units,'m/yr'),
+				dzostoga=dzostoga/1000;
+			end
 			time=self.time{i};time=time(1:end-1);
-			array=[dzostoga;time];
+			array=[dzostoga';time'];
 			arrays{i}=array;
 		end