Index: /issm/branches/trunk-larour-SLPS2020/src/m/contrib/larour/isintriangle.m
===================================================================
--- /issm/branches/trunk-larour-SLPS2020/src/m/contrib/larour/isintriangle.m	(revision 25775)
+++ /issm/branches/trunk-larour-SLPS2020/src/m/contrib/larour/isintriangle.m	(revision 25776)
@@ -10,4 +10,8 @@
 
 	a = y23 * dx + x32 * dy;
+	if isnan(a),
+		isin=0;
+		return;
+	end
 	if (a < minD || a > maxD)
 		isin=0; 
Index: /issm/branches/trunk-larour-SLPS2020/src/m/contrib/larour/ismip6.m
===================================================================
--- /issm/branches/trunk-larour-SLPS2020/src/m/contrib/larour/ismip6.m	(revision 25775)
+++ /issm/branches/trunk-larour-SLPS2020/src/m/contrib/larour/ismip6.m	(revision 25776)
@@ -246,6 +246,5 @@
 					for j=1:size(hg,2),
 						hgj=hg(:,j);
-						pos=find(hgj>.99); %we want fully grounded
-						%pos=find(hgj>0); %we want slightly grounded
+						pos=find(hgj>0); %we want slightly grounded
 						hgv(md.mesh.elements(pos,:),j)=1;
 					end
Index: /issm/branches/trunk-larour-SLPS2020/src/m/contrib/larour/oceanmip.m
===================================================================
--- /issm/branches/trunk-larour-SLPS2020/src/m/contrib/larour/oceanmip.m	(revision 25775)
+++ /issm/branches/trunk-larour-SLPS2020/src/m/contrib/larour/oceanmip.m	(revision 25776)
@@ -14,7 +14,9 @@
 		zostoga = {}; % global-mean thermosteric sea level anomalies (in mm) (time)
 		time = {}; %time in year
+		timem = {}; %time in months 
 		pbo = {}; %local ocean-bottom pressure changes. Also have zero mean over the oceans (lat,lon,time)
 		lat = {};
 		long = {};
+		mask = {};
 		scenario = {}; 
 		model = {}; 
@@ -44,7 +46,7 @@
 					%figure out time interval and remove historical:
 					time=ncread(file,'time');
+
 					pos=find(diff(time)<0); 
 					if isempty(pos),
-						%pos=(length(time)-12*100+1):length(time);
 						pos=1:length(time);
 					else
@@ -55,4 +57,5 @@
 					pos=pos(pos2);
 					time=time(pos); nt=length(time);
+					self.timem{end+1}=time;
 
 					%reduce datasets: 
@@ -67,5 +70,5 @@
 						zosm(:,:,year)=mean(zos(:,:,(i-11):i),3);
 					end
-					self.zos{end+1}=zosm; clear zos;
+					self.zos{end+1}=zosm; clear zos;  clear zosm;
 
 					%zostoga:
@@ -92,4 +95,7 @@
 					self.long{end+1}=ncread(file,'lon');
 
+					%mask: 
+					self.mask{end+1}=ncread(file,'slm'); 
+
 					%scenario and model: 
 					file=self.files{f};
@@ -102,5 +108,5 @@
 		end
 		%}}}
-		function self = interpolate(self,md) % {{{
+		function self = interpolate(self,md,slm) % {{{
 
 			%retrieve long and lat from mesh: 
@@ -110,4 +116,5 @@
 
 			for i=1:self.n,
+			%for i=1,
 				disp(['interpolating model ' self.model{i} ' onto model mesh']);
 
@@ -118,4 +125,5 @@
 					self.long{i}=self.long{i}(2:181,:,:);
 					self.lat{i}=self.lat{i}(2:181,:,:);
+					self.mask{i}=self.mask{i}(2:181,:,:);
 				end
 
@@ -125,6 +133,12 @@
 				pos=find(long<0); long(pos)=long(pos)+360;
 				lat=double(self.lat{i}); lat=lat(:);
-				[newl,uniqpos]=unique([lat,long],'rows','stable');
-				long=long(uniqpos); lat=lat(uniqpos);
+				if slm==1,
+					mask=double(self.mask{i}); mask=mask(:);
+					maskpos=find(mask==1);
+				else
+					maskpos=(1:length(long))';
+				end
+				[newl,uniqpos]=unique([lat(maskpos),long(maskpos)],'rows','stable');
+				long=long(maskpos(uniqpos)); lat=lat(maskpos(uniqpos));
 				index=delaunay(long,lat);
 
@@ -155,6 +169,6 @@
 					end
 
-					pbo=omip_pbo(:,:,j); pbo=pbo(:); pbo=pbo(uniqpos); 
-					zos=omip_zos(:,:,j); zos=zos(:); zos=zos(uniqpos);
+					pbo=omip_pbo(:,:,j); pbo=pbo(:); pbo=pbo(maskpos(uniqpos));
+					zos=omip_zos(:,:,j); zos=zos(:); zos=zos(maskpos(uniqpos));
 
 					pboj= InterpFromMeshToMesh2d(index,long,lat,pbo(newpos),meshlong,meshlat);