Changeset 13449 for issm/trunk-jpl/src/m/parameterization/setmask.py
- Timestamp:
- 09/26/12 12:06:42 (12 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
issm/trunk-jpl/src/m/parameterization/setmask.py
r13006 r13449 1 from numpy import * 2 import FlagElements as fe 1 import numpy 2 from FlagElements import * 3 3 4 4 def setmask(md, floatingicename, groundedicename): 5 #SETMASK - establish boundaries between grounded and floating ice. 6 # 7 # By default, ice is considered grounded. The contour floatingicename defines nodes 8 # for which ice is floating. The contour groundedicename defines nodes inside an floatingice, 9 # that are grounded (ie: ice rises, islands, etc ...) 10 # All input files are in the Argus format (extension .exp). 11 # 12 # Usage: 13 # md=setmask(md,floatingicename,groundedicename) 14 # 15 # Examples: 16 # md=setmask(md,'all',''); 17 # md=setmask(md,'Iceshelves.exp','Islands.exp'); 5 """ 6 SETMASK - establish boundaries between grounded and floating ice. 18 7 19 #%Get assigned fields 8 By default, ice is considered grounded. The contour floatingicename defines nodes 9 for which ice is floating. The contour groundedicename defines nodes inside an floatingice, 10 that are grounded (ie: ice rises, islands, etc ...) 11 All input files are in the Argus format (extension .exp). 12 13 Usage: 14 md=setmask(md,floatingicename,groundedicename) 15 16 Examples: 17 md=setmask(md,'all',''); 18 md=setmask(md,'Iceshelves.exp','Islands.exp'); 19 """ 20 21 #Get assigned fields 20 22 x = md.mesh.x 21 23 y = md.mesh.y … … 23 25 24 26 #Assign elementonfloatingice, elementongroundedice, vertexongroundedice and vertexonfloatingice. Only change at your own peril! This is synchronized heavily with the GroundingLineMigration module. {{{ 25 elementonfloatingice = fe.FlagElements(md, floatingicename)26 elementongroundedice = fe.FlagElements(md, groundedicename)27 elementonfloatingice = FlagElements(md, floatingicename) 28 elementongroundedice = FlagElements(md, groundedicename) 27 29 28 30 #Because groundedice nodes and elements can be included into an floatingice, we need to update. Remember, all the previous 29 31 #arrays come from domain outlines that can intersect one another: 30 32 31 elementonfloatingice = logical_and(elementonfloatingice,~elementongroundedice)32 elementongroundedice = ~elementonfloatingice33 elementonfloatingice = numpy.logical_and(elementonfloatingice,numpy.logical_not(elementongroundedice)) 34 elementongroundedice = numpy.logical_not(elementonfloatingice) 33 35 34 36 #the order here is important. we choose vertexongroundedice as default on the grounding line. 35 vertexonfloatingice = zeros(md.mesh.numberofvertices,'bool') 36 vertexongroundedice = zeros(md.mesh.numberofvertices,'bool') 37 38 pos=argwhere(elementongroundedice==1) 39 pos=md.mesh.elements[pos,:]-1 40 if pos.size: 41 vertexongroundedice[pos]=True 42 43 pos=argwhere(~vertexongroundedice) 44 if pos.size: 45 vertexonfloatingice[pos]=True; 46 #%}}} 37 vertexonfloatingice = numpy.zeros(md.mesh.numberofvertices,'bool') 38 vertexongroundedice = numpy.zeros(md.mesh.numberofvertices,'bool') 39 vertexongroundedice[md.mesh.elements[numpy.nonzero(elementongroundedice),:].astype(int)-1]=True 40 vertexonfloatingice[numpy.nonzero(numpy.logical_not(vertexongroundedice))]=True 41 #}}} 47 42 48 43 #Return: 49 md.mask.elementonfloatingice = double(elementonfloatingice) 50 md.mask.vertexonfloatingice = double(vertexonfloatingice) 51 md.mask.elementongroundedice = double(elementongroundedice) 52 md.mask.vertexongroundedice = double(vertexongroundedice) 53 md.mask.vertexonwater = zeros(md.mesh.numberofvertices) 54 md.mask.elementonwater = zeros(md.mesh.numberofelements) 44 md.mask.elementonfloatingice = elementonfloatingice.astype(float) 45 md.mask.vertexonfloatingice = vertexonfloatingice.astype(float) 46 md.mask.elementongroundedice = elementongroundedice.astype(float) 47 md.mask.vertexongroundedice = vertexongroundedice.astype(float) 48 md.mask.vertexonwater = numpy.zeros(md.mesh.numberofvertices) 49 md.mask.elementonwater = numpy.zeros(md.mesh.numberofelements) 50 55 51 return md
Note:
See TracChangeset
for help on using the changeset viewer.