0001 function md=geography(md,iceshelfname,icesheetname)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017 if ((nargin~=3) | (nargout~=1)),
0018 geographyusage();
0019 error('geography error message');
0020 end
0021 if ((~ischar(iceshelfname)) | (~ischar(icesheetname)) )
0022 geographyusage();
0023 error('geography error message');
0024 end
0025
0026 if md.counter>=2,
0027 choice=input('This model already has a geometry. Are you sure you want to go ahead? (y/n)','s');
0028 if ~strcmp(choice,'y')
0029 error('no geometry done ... exiting');
0030 end
0031 else
0032 if (md.counter~=1)
0033 error('geography error message: you need to run mesh.m first on this model');
0034 else
0035 md.counter=2;
0036 end
0037 end
0038
0039
0040 x=md.x;
0041 y=md.y;
0042 elements=md.elements;
0043
0044
0045 if strcmp(iceshelfname,''),
0046 gridoniceshelf=zeros(md.numberofgrids,1);
0047 elementoniceshelf=zeros(md.numberofelements,1);
0048 elseif strcmp(iceshelfname,'all'),
0049 gridoniceshelf=ones(md.numberofgrids,1);
0050 elementoniceshelf=ones(md.numberofelements,1);
0051 else
0052 [gridoniceshelf,elementoniceshelf]=ArgusContourToMesh(elements,x,y,expread(iceshelfname,1),'element and node',2);
0053 end
0054
0055 if strcmp(icesheetname,''),
0056 gridonicesheet=zeros(md.numberofgrids,1);
0057 elementonicesheet=zeros(md.numberofelements,1);
0058 else
0059 [gridonicesheet,elementonicesheet]=ArgusContourToMesh(elements,x,y,expread(icesheetname,1),'element and node',2);
0060 end
0061
0062
0063
0064 gridoniceshelf=double((gridoniceshelf & ~gridonicesheet));
0065 elementoniceshelf=double((elementoniceshelf & ~elementonicesheet));
0066 gridonicesheet=double(~gridoniceshelf);
0067 elementonicesheet=double(~elementoniceshelf);
0068
0069
0070 md.elementoniceshelf=elementoniceshelf;
0071 md.gridoniceshelf=gridoniceshelf;
0072
0073 md.elementonicesheet=elementonicesheet;
0074 md.gridonicesheet=gridonicesheet;
0075
0076
0077 if strcmp(iceshelfname,''), md.iceshelfoutline=''; elseif strcmp(iceshelfname,'all'), md.iceshelfoutline='all'; else md.iceshelfoutline=readfile(iceshelfname); end
0078 if strcmp(icesheetname,''), md.icesheetoutline=''; elseif strcmp(icesheetname,'all'), md.icesheetoutline='all'; else md.icesheetoutline=readfile(icesheetname); end
0079
0080 function geographyusage(),
0081 disp('INPUT md=geography(md,iceshelfname,icesheetname)');
0082 disp('geography: establish boundaries between grounded and floating ice.');
0083 disp(' By default, ice is considered grounded. The contour iceshelfname defines grids ');
0084 disp(' for which ice is floting. The contour icesheetname defines grids inside an iceshelf, ');
0085 disp(' that are grounded (ie: ice rises, islands, etc ...)');
0086 disp(' All input files are in the Argus format (extension .exp');