Changeset 9087

Timestamp:

07/20/11 15:00:12 (14 years ago)

Author:

schlegel

Message:

Add printflag

Location:

issm/trunk/test/NightlyRun

Files:

• test1501.m (modified) (2 diffs)
• test1502.m (modified) (2 diffs)

issm/trunk/test/NightlyRun/test1501.m

@@ +1 @@
+printingflag = false;
+
 md=mesh(model,'../Exp/Square.exp',350000);
 md=geography(md,'all','');
@@ -86 +88 @@
         PatchToVec(md.results.TransientSolution(2000).SurfaceMassBalance),...
         };
+
+if printingflag,
+
+        starttime = 360;
+        endtime = 2000;
+        res = 40;
+        ts = [starttime:res:endtime];
+
+        index = md.elements;
+        x1=md.x(index(:,1)); x2=md.x(index(:,2)); x3=md.x(index(:,3));
+        y1=md.y(index(:,1)); y2=md.y(index(:,2)); y3=md.y(index(:,3));
+        areas=(0.5*((x2-x1).*(y3-y1)-(y2-y1).*(x3-x1)));
+
+        thickness = [];
+        volume = [];
+        massbal = [];
+        velocity = [];
+        for t=starttime:endtime
+                thickness = [thickness PatchToVec(md.results.TransientSolution(t).Thickness)];
+                volume = [volume mean(md.results.TransientSolution(t).Thickness.value,2).*areas];
+                massbal = [massbal PatchToVec(md.results.TransientSolution(t).SurfaceMassBalance)];
+                velocity = [velocity PatchToVec(md.results.TransientSolution(t).Vel)];
+        end
+
+        figure('Position', [0 0 860 932])
+
+        options = plotoptions('data','transient_movie','unit','km');
+        options = options.list{1};
+        options = checkplotoptions(md,options);
+
+        %loop over the time steps
+        results=md.results.TransientSolution;
+        count = 1;
+        for i=ts
+
+                subplot(5,9,[28:31 37:40])
+                set(gca,'pos',get(gca,'pos')+[-0.08 -0.08 0.07 0.08])
+                field = 'Thickness';
+
+                %process data
+                [x y z elements is2d isplanet]=processmesh(md,results(i).(field),options);
+                [data datatype]=processdata(md,results(i).(field),options);
+
+                titlestring=[field ' at time ' num2str(results(i).time/md.yts) ' year'];
+                plot_unit(x,y,z,elements,data,is2d,isplanet,datatype,options)
+                options=changefieldvalue(options,'title',titlestring);
+                options=addfielddefault(options,'colorbar',1);
+                options=changefieldvalue(options,'caxis',[0 max(max(thickness))]);
+                applyoptions(md,[],options);
+
+                subplot(5,9,[33:36 42:45])
+                set(gca,'pos',get(gca,'pos')+[-0.00 -0.08 0.07 0.08])
+                field = 'Vel';
+
+                %process data
+                [x y z elements is2d isplanet]=processmesh(md,results(i).(field),options);
+                [data datatype]=processdata(md,results(i).(field),options);
+
+                titlestring=[field ' at time ' num2str(results(i).time/md.yts) ' year'];
+                plot_unit(x,y,z,elements,data,is2d,isplanet,datatype,options)
+                options=changefieldvalue(options,'title',titlestring);
+                options=addfielddefault(options,'colorbar',1);
+                options=changefieldvalue(options,'caxis',[0 max(max(velocity))]);
+                applyoptions(md,[],options);
+
+                subplot(5,4,1:4)
+                cla
+                set(gca,'pos',get(gca,'pos')+[-0.07 0.03 0.12 0.015])
+                plot(starttime:endtime,mean(massbal),'k','LineWidth', 4)
+                hold on
+                ya = ylim;
+                plot([i i], ya, 'r', 'LineWidth',6)
+                ylim(ya); xlim([starttime endtime]);
+                title('Surface Mass Balance','FontSize',14)
+                ylabel('m/year','FontSize',14)
+
+                subplot(5,4,5:8)
+                cla
+                set(gca,'pos',get(gca,'pos')+[-0.07 0.015 0.12 0.015])
+                plot(starttime:endtime,sum(volume)/1000/1000/1000,'LineWidth',4)
+                hold on
+                ya = ylim;
+                plot([i i], ya, 'r', 'LineWidth',6)
+                ylim(ya); xlim([starttime endtime]);
+                title('Ice Volume','FontSize',14)
+                ylabel('km^3','FontSize',14)
+
+                subplot(5,4,9:12)
+                cla
+                set(gca,'pos',get(gca,'pos')+[-0.07 0 0.12 0.015])
+                plot(starttime:endtime,mean(velocity)/1000, 'LineWidth', 4)
+                hold on
+                ya = ylim;
+                plot([i i], ya, 'r', 'LineWidth',6)
+                ylim(ya); xlim([starttime endtime]);
+                title('Mean Velocity','FontSize', 14)
+                ylabel('km/year','FontSize', 14)
+                xlabel('year','FontSize', 14)
+
+                set(gcf,'Renderer','zbuffer','color','white'); %fixes a bug on Mac OS X (not needed in future Matlab version)
+                if i==starttime,
+                        %initialize images and frame
+                        frame=getframe(gcf);
+                        [images,map]=rgb2ind(frame.cdata,256,'nodither');
+                        images(1,1,1,length(ts))=0;
+                else
+                        frame=getframe(gcf);
+                        images(:,:,1,count) = rgb2ind(frame.cdata,map,'nodither');
+                end
+
+                count = count+1;
+
+        end
+
+        filename='transawtooth2d.gif';
+        imwrite(images,map,filename,'DelayTime',1.0,'LoopCount',inf)
+
+end %printingflag

issm/trunk/test/NightlyRun/test1502.m

@@ +1 @@
+printingflag = false;
+
 md=mesh(model,'../Exp/Square.exp',700000);
 md=geography(md,'all','');
@@ -92 +94 @@
         PatchToVec(md.results.TransientSolution(2000).SurfaceMassBalance),...
         };
+
+if printingflag,
+
+        starttime = 360;
+        endtime = 2000;
+        res = 40;
+        ts = [starttime:res:endtime];
+
+        index = md.elements;
+        x1=md.x(index(:,1)); x2=md.x(index(:,2)); x3=md.x(index(:,3));
+        y1=md.y(index(:,1)); y2=md.y(index(:,2)); y3=md.y(index(:,3));
+        areas=(0.5*((x2-x1).*(y3-y1)-(y2-y1).*(x3-x1)));
+
+        thickness = [];
+        volume = [];
+        massbal = [];
+        velocity = [];
+        for t=starttime:endtime
+                thickness = [thickness PatchToVec(md.results.TransientSolution(t).Thickness)];
+                volume = [volume mean(md.results.TransientSolution(t).Thickness.value,2).*areas];
+                massbal = [massbal PatchToVec(md.results.TransientSolution(t).SurfaceMassBalance)];
+                velocity = [velocity PatchToVec(md.results.TransientSolution(t).Vel)];
+        end
+
+        figure('Position', [0 0 1060 1060])
+
+        options = plotoptions('data','transient_movie','unit','km');
+        options = options.list{1};
+        options = checkplotoptions(md,options);
+
+        %loop over the time steps
+        results=md.results.TransientSolution;
+        count = 1;
+        for i=ts
+
+                subplot(5,9,[28:31 37:40])
+                set(gca,'pos',get(gca,'pos')+[-0.08 -0.08 0.07 0.08])
+                field = 'Thickness';
+
+                %process data
+                [x y z elements is2d isplanet]=processmesh(md,results(i).(field),options);
+                [data datatype]=processdata(md,results(i).(field),options);
+
+                titlestring=[field ' at time ' num2str(results(i).time/md.yts) ' year'];
+                plot_unit(x,y,z,elements,data,is2d,isplanet,datatype,options)
+                options=changefieldvalue(options,'title',titlestring);
+                options=addfielddefault(options,'colorbar',1);
+                options=changefieldvalue(options,'caxis',[0 max(max(thickness))]);
+                applyoptions(md,[],options);
+
+                subplot(5,9,[33:36 42:45])
+                set(gca,'pos',get(gca,'pos')+[-0.01 -0.08 0.07 0.08])
+                field = 'Vel';
+
+                %process data
+                [x y z elements is2d isplanet]=processmesh(md,results(i).(field),options);
+                [data datatype]=processdata(md,results(i).(field),options);
+
+                titlestring=[field ' at time ' num2str(results(i).time/md.yts) ' year'];
+                plot_unit(x,y,z,elements,data,is2d,isplanet,datatype,options)
+                options=changefieldvalue(options,'title',titlestring);
+                options=addfielddefault(options,'colorbar',1);
+                options=changefieldvalue(options,'caxis',[0 max(max(velocity))]);
+                applyoptions(md,[],options);
+
+                subplot(5,4,1:4)
+                cla
+                set(gca,'pos',get(gca,'pos')+[-0.07 0.03 0.12 0.015])
+                plot(starttime:endtime,mean(massbal),'k','LineWidth', 4)
+                hold on
+                ya = ylim;
+                plot([i i], ya, 'r', 'LineWidth',6)
+                ylim(ya); xlim([starttime endtime]);
+                title('Surface Mass Balance','FontSize',14)
+                ylabel('m/year','FontSize',14)
+
+                subplot(5,4,5:8)
+                cla
+                set(gca,'pos',get(gca,'pos')+[-0.07 0.015 0.12 0.015])
+                plot(starttime:endtime,sum(volume)/1000/1000/1000,'LineWidth',4)
+                hold on
+                ya = ylim;
+                plot([i i], ya, 'r', 'LineWidth',6)
+                ylim(ya); xlim([starttime endtime]);
+                title('Ice Volume','FontSize',14)
+                ylabel('km^3','FontSize',14)
+
+                subplot(5,4,9:12)
+                cla
+                set(gca,'pos',get(gca,'pos')+[-0.07 0 0.12 0.015])
+                plot(starttime:endtime,mean(velocity)/1000, 'LineWidth', 4)
+                hold on
+                ya = ylim;
+                plot([i i], ya, 'r', 'LineWidth',6)
+                ylim(ya); xlim([starttime endtime]);
+                title('Mean Velocity','FontSize', 14)
+                ylabel('km/year','FontSize', 14)
+                xlabel('year','FontSize', 14)
+
+                set(gcf,'Renderer','zbuffer','color','white'); %fixes a bug on Mac OS X (not needed in future Matlab version)
+                if i==starttime,
+                        %initialize images and frame
+                        frame=getframe(gcf);
+                        [images,map]=rgb2ind(frame.cdata,256,'nodither');
+                        images(1,1,1,length(ts))=0;
+                else
+                        frame=getframe(gcf);
+                        images(:,:,1,count) = rgb2ind(frame.cdata,map,'nodither');
+                end
+
+                count = count+1;
+
+        end
+
+        filename='transawtooth3d.gif';
+        imwrite(images,map,filename,'DelayTime',1.0,'LoopCount',inf)
+
+end %printingflag