Index: /issm/trunk-jpl/src/m/classes/mesh3dprisms.m
===================================================================
--- /issm/trunk-jpl/src/m/classes/mesh3dprisms.m	(revision 21416)
+++ /issm/trunk-jpl/src/m/classes/mesh3dprisms.m	(revision 21417)
@@ -104,7 +104,18 @@
 			md = checkfield(md,'fieldname','mesh.vertexonbase','size',[md.mesh.numberofvertices 1],'values',[0 1]);
 			md = checkfield(md,'fieldname','mesh.vertexonsurface','size',[md.mesh.numberofvertices 1],'values',[0 1]);
+			md = checkfield(md,'fieldname','mesh.average_vertex_connectivity','>=',24,'message','''mesh.average_vertex_connectivity'' should be at least 24 in 3d');
+
+			%Check that mesh follows the geometry
 			md = checkfield(md,'fieldname','mesh.z','>=',md.geometry.base-10^-10,'message','''mesh.z'' lower than bedrock');
 			md = checkfield(md,'fieldname','mesh.z','<=',md.geometry.surface+10^-10,'message','''mesh.z'' higher than surface elevation');
-			md = checkfield(md,'fieldname','mesh.average_vertex_connectivity','>=',24,'message','''mesh.average_vertex_connectivity'' should be at least 24 in 3d');
+			if any(max(abs(project2d(md,md.mesh.z,1)-project2d(md,md.geometry.base,1)))>1e-11),
+				md = checkmessage(md,'md.mesh.z is not consistent with md.geometry.base, you changed the geometry after extrusion');
+			end
+			if any(max(abs(project2d(md,md.mesh.z,md.mesh.numberoflayers)-project2d(md,md.geometry.surface,1)))>1e-11),
+				md = checkmessage(md,'md.mesh.z is not consistent with md.geometry.surface, you changed the geometry after extrusion !!');
+			end
+			if any(max(abs(project2d(md,md.mesh.z,md.mesh.numberoflayers)-project2d(md,md.mesh.z,1) - project2d(md,md.geometry.thickness,1)))>1e-11),
+				md = checkmessage(md,'md.mesh.z is not consistent with md.geometry.thickness, you changed the geometry after extrusion !!');
+			end
 		end % }}}
 		function disp(self) % {{{
Index: /issm/trunk-jpl/src/m/classes/model.m
===================================================================
--- /issm/trunk-jpl/src/m/classes/model.m	(revision 21416)
+++ /issm/trunk-jpl/src/m/classes/model.m	(revision 21417)
@@ -223,4 +223,6 @@
 			if ~isnan(md.initialization.epl_head),md.initialization.epl_head=project2d(md,md.initialization.epl_head,1);end;
 			if ~isnan(md.initialization.epl_thickness),md.initialization.epl_thickness=project2d(md,md.initialization.epl_thickness,1);end;
+			if ~isnan(md.initialization.waterfraction),md.initialization.waterfraction=project2d(md,md.initialization.waterfraction,1);end;
+			if ~isnan(md.initialization.watercolumn),md.initialization.watercolumn=project2d(md,md.initialization.watercolumn,1);end;
 			%gia
 			if ~isnan(md.gia.mantle_viscosity), md.gia.mantle_viscosity=project2d(md,md.gia.mantle_viscosity,1); end
@@ -244,4 +246,5 @@
 			md.masstransport.spcthickness=project2d(md,md.masstransport.spcthickness,md.mesh.numberoflayers);
 			if ~isnan(md.damage.spcdamage), md.damage.spcdamage=project2d(md,md.damage.spcdamage,md.mesh.numberoflayers); end
+			if numel(md.levelset.spclevelset)>1, md.levelset.spclevelset=project2d(md,md.levelset.spclevelset,md.mesh.numberoflayers); end
 			md.thermal.spctemperature=project2d(md,md.thermal.spctemperature,md.mesh.numberoflayers);
 
@@ -270,8 +273,15 @@
 				md.basalforcings.groundedice_melting_rate=project2d(md,md.basalforcings.groundedice_melting_rate,1); 
 			end
-			if ~isnan(md.basalforcings.floatingice_melting_rate),
+			if isprop(md.basalforcings,'floatingice_melting_rate') & ~isnan(md.basalforcings.floatingice_melting_rate),
 				md.basalforcings.floatingice_melting_rate=project2d(md,md.basalforcings.floatingice_melting_rate,1); 
 			end
 			md.basalforcings.geothermalflux=project2d(md,md.basalforcings.geothermalflux,1); %bedrock only gets geothermal flux
+
+			if isprop(md.calving,'coeff') & ~isnan(md.calving.coeff),
+				md.calving.coeff=project2d(md,md.calving.coeff,1); 
+			end
+			if isprop(md.calving,'meltingrate') & ~isnan(md.calving.meltingrate),
+				md.calving.meltingrate=project2d(md,md.calving.meltingrate,1); 
+			end
 
 			%update of connectivity matrix