[12888] | 1 | import numpy
|
---|
| 2 | import os
|
---|
| 3 | #from basinzoom import *
|
---|
| 4 | #from ContourToMesh import *
|
---|
| 5 | from MatlabFuncs import *
|
---|
| 6 |
|
---|
[12112] | 7 | def FlagElements(md,region):
|
---|
[12888] | 8 | """
|
---|
| 9 | FLAGELEMENTS - flag the elements in an region
|
---|
[12112] | 10 |
|
---|
[12944] | 11 | The region can be given with an exp file, a list of elements.
|
---|
[12888] | 12 |
|
---|
[12944] | 13 | Usage:
|
---|
| 14 | flag=FlagElements(md,region);
|
---|
[12888] | 15 |
|
---|
[12944] | 16 | Example:
|
---|
| 17 | flag=FlagElements(md,'all');
|
---|
| 18 | flag=FlagElements(md,'');
|
---|
| 19 | flag=FlagElements(md,'Domain.exp');
|
---|
| 20 | flag=FlagElements(md,'~Domain.exp');
|
---|
| 21 | flag=FlagElements(md,md.mask.elementongroundedice);
|
---|
[12888] | 22 | """
|
---|
| 23 |
|
---|
[13098] | 24 | if isinstance(region,(str,unicode)):
|
---|
[12888] | 25 | if not region:
|
---|
| 26 | flag=numpy.zeros(md.mesh.numberofelements,'bool')
|
---|
| 27 | invert=0
|
---|
| 28 | elif strcmpi(region,'all'):
|
---|
| 29 | flag=numpy.ones(md.mesh.numberofelements,'bool')
|
---|
| 30 | invert=0
|
---|
[12112] | 31 | else:
|
---|
| 32 | #make sure that we actually don't want the elements outside the domain outline!
|
---|
[12888] | 33 | if strcmpi(region[0],'~'):
|
---|
[12112] | 34 | region=region[1:]
|
---|
[12888] | 35 | invert=1
|
---|
[12112] | 36 | else:
|
---|
[12888] | 37 | invert=0
|
---|
| 38 |
|
---|
[12112] | 39 | #does the region domain outline exist or do we have to look for xlim,ylim in basinzoom?
|
---|
[12888] | 40 | if not os.path.exists(region):
|
---|
| 41 | if len(region)>3 and not strcmp(region[-4:],'.exp'):
|
---|
| 42 | raise IOError("Error: File 'region' not found!" % region)
|
---|
| 43 | xlim,ylim=basinzoom('basin',region)
|
---|
| 44 | flag_nodes=numpy.logical_and(numpy.logical_and(md.mesh.x<xlim[1],md.mesh.x>xlim[0]),numpy.logical_and(md.mesh.y<ylim[1],md.mesh.y>ylim[0])).astype(float)
|
---|
| 45 | flag=numpy.prod(flag_nodes[md.mesh.elements],axis=1)
|
---|
[12112] | 46 | else:
|
---|
| 47 | #ok, flag elements
|
---|
[12888] | 48 | flag=ContourToMesh(md.mesh.elements[:,0:3],md.mesh.x,md.mesh.y,region,'element',1)
|
---|
| 49 |
|
---|
[12112] | 50 | if invert:
|
---|
[12888] | 51 | flag=numpy.logical_not(flag)
|
---|
| 52 |
|
---|
| 53 | elif isinstance(region,numpy.nparray) or isinstance(region,bool):
|
---|
| 54 | if not numpy.size(region,0)==md.mesh.numberofelements:
|
---|
| 55 | raise TypeError("Flaglist for region must be of same size as number of elements in model.")
|
---|
| 56 | flag=region
|
---|
| 57 |
|
---|
[12112] | 58 | else:
|
---|
[12888] | 59 | raise TypeError("Invalid region option")
|
---|
[12112] | 60 |
|
---|
[12888] | 61 | return flag
|
---|
| 62 |
|
---|