source: issm/trunk/src/m/geometry/FlagElements.py@ 13975

import numpy
import os
#from basinzoom import *
from ContourToMesh import *
from MatlabFuncs import *

def FlagElements(md,region):
        """
        FLAGELEMENTS - flag the elements in an region

           The region can be given with an exp file, a list of elements.

           Usage: 
              flag=FlagElements(md,region);

           Example:
              flag=FlagElements(md,'all');
              flag=FlagElements(md,'');
              flag=FlagElements(md,'Domain.exp');
              flag=FlagElements(md,'~Domain.exp');
              flag=FlagElements(md,md.mask.elementongroundedice);
        """

        if   isinstance(region,(str,unicode)):
                if   not region:
                        flag=numpy.zeros(md.mesh.numberofelements,'bool')
                        invert=0
                elif strcmpi(region,'all'):
                        flag=numpy.ones(md.mesh.numberofelements,'bool')
                        invert=0
                else:
                        #make sure that we actually don't want the elements outside the domain outline!
                        if strcmpi(region[0],'~'):
                                region=region[1:]
                                invert=1
                        else:
                                invert=0

                        #does the region domain outline exist or do we have to look for xlim,ylim in basinzoom?
                        if not os.path.exists(region):
                                if len(region)>3 and not strcmp(region[-4:],'.exp'):
                                        raise IOError("Error: File 'region' not found!" % region)
                                raise RuntimeError("FlagElements.py calling basinzoom.py is not complete.")
                                xlim,ylim=basinzoom('basin',region)
                                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)
                                flag=numpy.prod(flag_nodes[md.mesh.elements],axis=1)
                        else:
                                #ok, flag elements
                                [flag,dum]=ContourToMesh(md.mesh.elements[:,0:3].copy(),md.mesh.x.reshape(-1,1),md.mesh.y.reshape(-1,1),region,'element',1)

                if invert:
                        flag=numpy.logical_not(flag)

        elif isinstance(region,numpy.ndarray) or isinstance(region,bool):
                if not numpy.size(region,0)==md.mesh.numberofelements:
                        raise TypeError("Flaglist for region must be of same size as number of elements in model.")
                flag=region

        else:
                raise TypeError("Invalid region option")

        return flag