Index: /issm/trunk-jpl/examples/SlrGRACE/runme.m
===================================================================
--- /issm/trunk-jpl/examples/SlrGRACE/runme.m	(revision 26258)
+++ /issm/trunk-jpl/examples/SlrGRACE/runme.m	(revision 26259)
@@ -2,5 +2,5 @@
 addpath('../Data','../Functions');
 
-steps=[1:7];
+steps=[6:7];
 
 if any(steps==1) % {{{ 
@@ -142,5 +142,5 @@
 	sol2 = (md.results.TransientSolution.Sealevel-md.results.TransientSolution.Bed)*1000;	% [mm]
 
-	sol_name={'Change in water equivalent height [cm]', 'Relative sea-level [mm]'};
+	sol_name={'Change in water equivalent height [cm]', 'Relative sea level [mm]'};
 	fig_name={'Fig_dH.pdf','Fig_slr.pdf'};
 
@@ -200,27 +200,34 @@
 	disp('Projecting  loads onto the mesh...');
 	time_range = 2007 + [15 350]/365;
-	water_load = grace(md.mesh.elements,md.mesh.lat,md.mesh.long,time_range(1),time_range(2));
-
-	[ne,nt]=size(water_load);
-	md.solidearth.surfaceload.icethicknesschange = zeros(ne+1,nt);
-	md.solidearth.surfaceload.icethicknesschange(1:ne,:) = water_load*md.materials.rho_freshwater/md.materials.rho_ice;
-	md.solidearth.surfaceload.icethicknesschange(ne+1,:)=[1:nt]; % times
-
+	onvertex = 1; % map data on vertex. If 0, it maps on the elemental centroid. 
+	water_load = grace(md.mesh.elements,md.mesh.lat,md.mesh.long,time_range(1),time_range(2),onvertex);
+	rho_water2ice = md.materials.rho_freshwater/md.materials.rho_ice; 
+	ice_load = water_load*rho_water2ice; % ice height equivalent. 
+
+	% masstransport evalulates diff between the successive times, so we should cumsum. 
+	num_time = size(ice_load,2); 
+	md.masstransport.spcthickness = [md.geometry.thickness + ice_load];
+	md.masstransport.spcthickness(md.mesh.numberofvertices+1,:)=[0:num_time-1]; 
+
+	%Physics 
 	md.transient.issmb=0;
-	md.transient.ismasstransport=0;
 	md.transient.isstressbalance=0;
 	md.transient.isthermal=0;
-	md.transient.isgia=1;	md.transient.isslr=1;
-
-	md.timestepping.start_time=-1;
-	md.timestepping.final_time=nt;
+	md.transient.ismasstransport=1;
+	md.transient.isslc=1;
+
+	md.solidearth.settings.rigid=1;
+	md.solidearth.settings.elastic=1;
+	md.solidearth.settings.rotation=1;
+	
+	%time stepping:
+	md.timestepping.start_time=0;
 	md.timestepping.time_step=1;
-	md.timestepping.interp_forcing=0;
-	md.settings.output_frequency=1;
+	md.timestepping.final_time=num_time; 
 
 	md.cluster=generic('name',oshostname(),'np',3);
 	md.verbose=verbose('111111111');
 
-	md.solidearth.requested_outputs = {'Sealevel','SealevelRSL'};
+	md.solidearth.requested_outputs = {'Sealevel','Bed'};
 
 	md=solve(md,'tr');
@@ -232,14 +239,14 @@
 	md = loadmodel('./Models/SlrGRACE_Transient');
 
-	time = md.solidearth.surfaceload.icethicknesschange(end,:);
+	time = md.masstransport.spcthickness(end,:);
 
 	for tt=1:length(time)
-		gmsl(tt) = md.results.TransientSolution(tt).Bslr*1000; % GMSL rate mm/yr
-		sol1(:,tt) = md.solidearth.surfaceload.icethicknesschange(1:end-1,tt)*100;	% ice equivalent height [cm/yr]
-		sol2(:,tt) = md.results.TransientSolution(tt+1).SealevelRSL*1000;	% mm/yr
+		gmsl(tt) = md.results.TransientSolution(tt).Bslc*1000; % [mm]
+		sol1(:,tt) = (md.masstransport.spcthickness(1:end-1,tt)-md.geometry.thickness)*100; % [cm]
+		sol2(:,tt) = (md.results.TransientSolution(tt).Sealevel-md.results.TransientSolution(tt).Bed)*1000;	% [mm]
 	end
-	sol_name = {'Change in water equivalent height [cm]', 'Relative sea-level [mm]'};
+	sol_name = {'Change in water equivalent height [cm]', 'Relative sea level [mm]'};
 	movie_name = {'Movie_dH.mp4','Movie_slr.mp4'};
-
+	
 	res = 1.0;
 
@@ -284,5 +291,5 @@
 				geoshow(flipud(coastlat),flipud(coastlon),'DisplayType','polygon','FaceColor','white');
 			else
-				geoshow(coastlat,coastlon,'DisplayType','polygon','FaceColor','white');
+				geoshow(coastlat,coastlon,'DisplayType','polygon','FaceColor','none');
 			end
 			plot(coastlon,coastlat,'k'); hold off;