def fields(depth, Rp, Rm, Gam, lrec, lsrc, zsrc, TM):
    r"""Calculate Pu+, Pu-, Pd+, Pd-.

    This is a modified version of empymod.kernel.fields(). See the original
    version for more information.

    """
    # Booleans if src in first or last layer; swapped if up=True
    first_layer = lsrc == 0
    last_layer = lsrc == depth.size-1

    # Depths; dp and dm are swapped if up=True
    if lsrc != depth.size-1:
        ds = depth[lsrc+1]-depth[lsrc]
        dp = depth[lsrc+1]-zsrc
    dm = zsrc-depth[lsrc]

    # Rm and Rp; swapped if up=True
    Rmp = Rm
    Rpm = Rp

    # Boolean if plus or minus has to be calculated
    if TM:
        plus = False
    else:
        plus = True

    # Sign-switches
    pm = 1     # + if plus=True, - if plus=False
    if not plus:
        pm = -1

    # Calculate down- and up-going fields
    for up in [False, True]:

        # No upgoing field if rec is in last layer or below src
        if up and (lrec == depth.size-1 or lrec > lsrc):
            Puu = np.full_like(Gam[:, :, lsrc, :], 0+0j)
            Pud = np.full_like(Gam[:, :, lsrc, :], 0+0j)
            continue
        # No downgoing field if rec is in first layer or above src
        if not up and (lrec == 0 or lrec < lsrc):
            Pdu = np.full_like(Gam[:, :, lsrc, :], 0+0j)
            Pdd = np.full_like(Gam[:, :, lsrc, :], 0+0j)
            continue

        # Swaps if up=True
        if up:
            dp, dm = dm, dp
            Rmp, Rpm = Rpm, Rmp
            first_layer, last_layer = last_layer, first_layer

        # Calculate Pu+, Pu-, Pd+, Pd-; rec in src layer; Eqs  81/82, A-8/A-9
        iGam = Gam[:, :, lsrc, :]
        if last_layer:  # If src/rec are in top (up) or bottom (down) layer
            Pd = Rmp*np.exp(-iGam*dm)
            Pu = np.full_like(Gam[:, :, lsrc, :], 0+0j)
        else:           # If src and rec are in any layer in between
            Ms = 1 - Rmp*Rpm*np.exp(-2*iGam*ds)
            Pd = Rmp/Ms*np.exp(-iGam*dm)
            Pu = Rmp/Ms*pm*Rpm*np.exp(-iGam*(ds+dp))

        # Store P's
        if up:
            Puu = Pu
            Pud = Pd
        else:
            Pdu = Pd
            Pdd = Pu

    # Return fields (up- and downgoing)
    return Puu, Pud, Pdu, Pdd