<commit_before>from .matutils import mkvc, ndgrid
import numpy as np

def surface2ind_topo(mesh, topo, gridLoc='CC'):
# def genActiveindfromTopo(mesh, topo):
    """
    Get active indices from topography
    """


    if mesh.dim == 3:
        from scipy.interpolate import NearestNDInterpolator
        Ftopo = NearestNDInterpolator(topo[:,:2], topo[:,2])

        if gridLoc == 'CC':
            XY = ndgrid(mesh.vectorCCx, mesh.vectorCCy)
            Zcc = mesh.gridCC[:,2].reshape((np.prod(mesh.vnC[:2]), mesh.nCz), order='F')

            gridTopo = Ftopo(XY)
            actind = [gridTopo[ixy] <= Zcc[ixy,:] for ixy in range(np.prod(mesh.vnC[0]))]
            actind = np.hstack(actind)

        elif gridLoc == 'N':

            XY = ndgrid(mesh.vectorNx, mesh.vectorNy)
            gridTopo = Ftopo(XY).reshape(mesh.vnN[:2], order='F')

            if mesh._meshType not in ['TENSOR', 'CYL', 'BASETENSOR']:
                raise NotImplementedError('Nodal surface2ind_topo not implemented for {0!s} mesh'.format(mesh._meshType))

            Nz = mesh.vectorNz[1:] # TODO: this will only work for tensor meshes
            actind = np.array([False]*mesh.nC).reshape(mesh.vnC, order='F')

            for ii in range(mesh.nCx):
                for jj in range(mesh.nCy):
                     actind[ii,jj,:] = [np.all(gridTopo[ii:ii+2, jj:jj+2] >= Nz[kk]) for kk in range(len(Nz)) ]

    elif mesh.dim == 2:
        from scipy.interpolate import interp1d
        Ftopo = interp1d(topo[:,0], topo[:,1])

        if gridLoc == 'CC':
            gridTopo = Ftopo(mesh.gridCC[:,0])
            actind = mesh.gridCC[:,1] <= gridTopo

        elif gridLoc == 'N':

            gridTopo = Ftopo(mesh.vectorNx)
            if mesh._meshType not in ['TENSOR', 'CYL', 'BASETENSOR']:
                raise NotImplementedError('Nodal surface2ind_topo not implemented for {0!s} mesh'.format(mesh._meshType))

            Ny = mesh.vectorNy[1:] # TODO: this will only work for tensor meshes
            actind = np.array([False]*mesh.nC).reshape(mesh.vnC, order='F')

            for ii in range(mesh.nCx):
                actind[ii,:] = [np.all(gridTopo[ii:ii+2] > Ny[kk]) for kk in range(len(Ny)) ]

    else:
        raise NotImplementedError('surface2ind_topo not implemented for 1D mesh')

    return mkvc(actind)

<commit_msg>change interpolation flag to "extrapolate" for surf2indtopo<commit_after>from .matutils import mkvc, ndgrid
import numpy as np

def surface2ind_topo(mesh, topo, gridLoc='CC'):
# def genActiveindfromTopo(mesh, topo):
    """
    Get active indices from topography
    """


    if mesh.dim == 3:
        from scipy.interpolate import NearestNDInterpolator
        Ftopo = NearestNDInterpolator(topo[:,:2], topo[:,2])

        if gridLoc == 'CC':
            XY = ndgrid(mesh.vectorCCx, mesh.vectorCCy)
            Zcc = mesh.gridCC[:,2].reshape((np.prod(mesh.vnC[:2]), mesh.nCz), order='F')

            gridTopo = Ftopo(XY)
            actind = [gridTopo[ixy] <= Zcc[ixy,:] for ixy in range(np.prod(mesh.vnC[0]))]
            actind = np.hstack(actind)

        elif gridLoc == 'N':

            XY = ndgrid(mesh.vectorNx, mesh.vectorNy)
            gridTopo = Ftopo(XY).reshape(mesh.vnN[:2], order='F')

            if mesh._meshType not in ['TENSOR', 'CYL', 'BASETENSOR']:
                raise NotImplementedError('Nodal surface2ind_topo not implemented for {0!s} mesh'.format(mesh._meshType))

            Nz = mesh.vectorNz[1:] # TODO: this will only work for tensor meshes
            actind = np.array([False]*mesh.nC).reshape(mesh.vnC, order='F')

            for ii in range(mesh.nCx):
                for jj in range(mesh.nCy):
                     actind[ii,jj,:] = [np.all(gridTopo[ii:ii+2, jj:jj+2] >= Nz[kk]) for kk in range(len(Nz)) ]

    elif mesh.dim == 2:
        from scipy.interpolate import interp1d
        Ftopo = interp1d(topo[:,0], topo[:,1], fill_value="extrapolate")

        if gridLoc == 'CC':
            gridTopo = Ftopo(mesh.gridCC[:,0])
            actind = mesh.gridCC[:,1] <= gridTopo

        elif gridLoc == 'N':

            gridTopo = Ftopo(mesh.vectorNx)
            if mesh._meshType not in ['TENSOR', 'CYL', 'BASETENSOR']:
                raise NotImplementedError('Nodal surface2ind_topo not implemented for {0!s} mesh'.format(mesh._meshType))

            Ny = mesh.vectorNy[1:] # TODO: this will only work for tensor meshes
            actind = np.array([False]*mesh.nC).reshape(mesh.vnC, order='F')

            for ii in range(mesh.nCx):
                actind[ii,:] = [np.all(gridTopo[ii:ii+2] > Ny[kk]) for kk in range(len(Ny)) ]

    else:
        raise NotImplementedError('surface2ind_topo not implemented for 1D mesh')

    return mkvc(actind)


{'fraction_non_alphanumeric': 0.14219212865002115, 'fraction_numerical': 0.011214557765552264, 'mean_word_length': 2.4808541973490428, 'pattern_counts': {'":': 0, '<': 7, '<?xml version=': 0, '>': 7, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '5293455', 'n_tokens_mistral': 1708, 'n_tokens_neox': 1635, 'n_words': 361}