def fem(ab, off, angle, zsrc, zrec, lsrc, lrec, depth, freq, etaH, etaV, zetaH,
        zetaV, xdirect, isfullspace, ht, htarg, use_ne_eval, msrc, mrec,
        loop_freq, loop_off, conv=True):
    r"""Return the electromagnetic frequency-domain response.

    This function is called from one of the above modelling routines. No
    input-check is carried out here. See the main description of :mod:`model`
    for information regarding input and output parameters.

    This function can be directly used if you are sure the provided input is in
    the correct format. This is useful for inversion routines and similar, as
    it can speed-up the calculation by omitting input-checks.

    """
    # Preallocate array
    fEM = np.zeros((freq.size, off.size), dtype=complex)

    # Initialize kernel count
    # (how many times the wavenumber-domain kernel was calld)
    kcount = 0

    # If <ab> = 36 (or 63), fEM-field is zero
    if ab in [36, ]:
        return fEM, kcount, conv

    # Get full-space-solution if xdirect=True and model is a full-space or
    # if src and rec are in the same layer.
    if xdirect and (isfullspace or lsrc == lrec):
        fEM += kernel.fullspace(off, angle, zsrc, zrec, etaH[:, lrec],
                                etaV[:, lrec], zetaH[:, lrec], zetaV[:, lrec],
                                ab, msrc, mrec)

    # If `xdirect = None` we set it here to True, so it is NOT calculated in
    # the wavenumber domain. (Only reflected fields are returned.)
    if xdirect is None:
        xdir = True
    else:
        xdir = xdirect

    # Get angle dependent factors
    factAng = kernel.angle_factor(angle, ab, msrc, mrec)

    # Compute required lambdas for given hankel-filter-base
    # This should be in utils, but this is a backwards-incompatible change.
    # Move this to utils for version 2.0.
    if ht == 'fht':
        # htarg[0] = filter; htarg[1] = pts_per_dec
        lambd, int_pts = transform.get_spline_values(htarg[0], off, htarg[1])
        if not loop_off:
            htarg = (htarg[0], htarg[1], lambd, int_pts)

    # If not full-space with xdirect calculate fEM-field
    if not isfullspace*xdir:
        calc = getattr(transform, ht)
        if loop_freq:

            for i in range(freq.size):
                out = calc(zsrc, zrec, lsrc, lrec, off, factAng, depth, ab,
                           etaH[None, i, :], etaV[None, i, :],
                           zetaH[None, i, :], zetaV[None, i, :], xdir,
                           htarg, use_ne_eval, msrc, mrec)
                fEM[None, i, :] += out[0]
                kcount += out[1]
                conv *= out[2]

        elif loop_off:
            for i in range(off.size):

                # See comments above where it says "ht == 'fht'".
                # Get pre-calculated lambd, int_pts for this offset
                if ht == 'fht':
                    htarg = (htarg[0], htarg[1], lambd[None, i, :], int_pts[i])

                out = calc(zsrc, zrec, lsrc, lrec, off[None, i],
                           factAng[None, i], depth, ab, etaH, etaV, zetaH,
                           zetaV, xdir, htarg, use_ne_eval, msrc, mrec)
                fEM[:, None, i] += out[0]
                kcount += out[1]
                conv *= out[2]
        else:
            out = calc(zsrc, zrec, lsrc, lrec, off, factAng, depth, ab, etaH,
                       etaV, zetaH, zetaV, xdir, htarg, use_ne_eval, msrc,
                       mrec)
            fEM += out[0]
            kcount += out[1]
            conv *= out[2]

    return fEM, kcount, conv
def fem ( ab , off , angle , zsrc , zrec , lsrc , lrec , depth , freq , etaH , etaV , zetaH , zetaV , xdirect , isfullspace , ht , htarg , use_ne_eval , msrc , mrec , loop_freq , loop_off , conv = True ) : fEM = np . zeros ( ( freq . size , off . size ) , dtype = complex ) kcount = 0 if ab in [ 36 , ] : return fEM , kcount , conv if xdirect and ( isfullspace or lsrc == lrec ) : fEM += kernel . fullspace ( off , angle , zsrc , zrec , etaH [ : , lrec ] , etaV [ : , lrec ] , zetaH [ : , lrec ] , zetaV [ : , lrec ] , ab , msrc , mrec ) if xdirect is None : xdir = True else : xdir = xdirect factAng = kernel . angle_factor ( angle , ab , msrc , mrec ) if ht == 'fht' : lambd , int_pts = transform . get_spline_values ( htarg [ 0 ] , off , htarg [ 1 ] ) if not loop_off : htarg = ( htarg [ 0 ] , htarg [ 1 ] , lambd , int_pts ) if not isfullspace * xdir : calc = getattr ( transform , ht ) if loop_freq : for i in range ( freq . size ) : out = calc ( zsrc , zrec , lsrc , lrec , off , factAng , depth , ab , etaH [ None , i , : ] , etaV [ None , i , : ] , zetaH [ None , i , : ] , zetaV [ None , i , : ] , xdir , htarg , use_ne_eval , msrc , mrec ) fEM [ None , i , : ] += out [ 0 ] kcount += out [ 1 ] conv *= out [ 2 ] elif loop_off : for i in range ( off . size ) : if ht == 'fht' : htarg = ( htarg [ 0 ] , htarg [ 1 ] , lambd [ None , i , : ] , int_pts [ i ] ) out = calc ( zsrc , zrec , lsrc , lrec , off [ None , i ] , factAng [ None , i ] , depth , ab , etaH , etaV , zetaH , zetaV , xdir , htarg , use_ne_eval , msrc , mrec ) fEM [ : , None , i ] += out [ 0 ] kcount += out [ 1 ] conv *= out [ 2 ] else : out = calc ( zsrc , zrec , lsrc , lrec , off , factAng , depth , ab , etaH , etaV , zetaH , zetaV , xdir , htarg , use_ne_eval , msrc , mrec ) fEM += out [ 0 ] kcount += out [ 1 ] conv *= out [ 2 ] return fEM , kcount , conv
r"""Return the electromagnetic frequency-domain response.

    This function is called from one of the above modelling routines. No
    input-check is carried out here. See the main description of :mod:`model`
    for information regarding input and output parameters.

    This function can be directly used if you are sure the provided input is in
    the correct format. This is useful for inversion routines and similar, as
    it can speed-up the calculation by omitting input-checks.