How to use the plasmapy.formulary.parameters.thermal_speed function in plasmapy

v = v.to(u.m / u.s)
        # Catching case where drift velocity has default value, and
        # needs to be assigned units
        v_drift = _v_drift_units(v_drift)
        # convert temperature to Kelvins
        T = T.to(u.K, equivalencies=u.temperature_energy())
        if np.isnan(vTh):
            # get thermal velocity and thermal velocity squared
            vTh = parameters.thermal_speed(T, particle=particle, method="most_probable")
        elif not np.isnan(vTh):
            # check units of thermal velocity
            vTh = vTh.to(u.m / u.s)
    elif np.isnan(vTh) and units == "unitless":
        # assuming unitless temperature is in Kelvins
        vTh = (
            parameters.thermal_speed(T * u.K, particle=particle, method="most_probable")
        ).si.value
    # getting square of thermal speed
    vThSq = vTh ** 2
    # get square of relative particle speed
    vSq = (v - v_drift) ** 2
    # calculating distribution function
    coeff1 = (np.pi * vThSq) ** (-1)
    coeff2 = 2 * np.pi * (v - v_drift)
    expTerm = np.exp(-vSq / vThSq)
    distFunc = coeff1 * coeff2 * expTerm
    if units == "units":
        return distFunc.to(u.s / u.m)
    elif units == "unitless":
        return distFunc

v = v.to(u.m / u.s)
        # Catching case where drift velocity has default value, and
        # needs to be assigned units
        v_drift = _v_drift_units(v_drift)
        # convert temperature to Kelvins
        T = T.to(u.K, equivalencies=u.temperature_energy())
        if np.isnan(vTh):
            # get thermal velocity and thermal velocity squared
            vTh = parameters.thermal_speed(T, particle=particle, method="most_probable")
        elif not np.isnan(vTh):
            # check units of thermal velocity
            vTh = vTh.to(u.m / u.s)
    elif np.isnan(vTh) and units == "unitless":
        # assuming unitless temperature is in Kelvins
        vTh = (
            parameters.thermal_speed(T * u.K, particle=particle, method="most_probable")
        ).si.value
    # getting square of thermal speed
    vThSq = vTh ** 2
    # get square of relative particle speed
    vSq = (v - v_drift) ** 2
    # calculating distribution function
    coeff1 = (np.pi * vThSq) ** (-3 / 2)
    coeff2 = 4 * np.pi * vSq
    expTerm = np.exp(-vSq / vThSq)
    distFunc = coeff1 * coeff2 * expTerm
    if units == "units":
        return distFunc.to(u.s / u.m)
    elif units == "unitless":
        return distFunc

v = v.to(u.m / u.s)
        # Catching case where drift velocities have default values, they
        # need to be assigned units
        v_drift = _v_drift_units(v_drift)
        # convert temperature to Kelvins
        T = T.to(u.K, equivalencies=u.temperature_energy())
        if np.isnan(vTh):
            # get thermal velocity and thermal velocity squared
            vTh = parameters.thermal_speed(T, particle=particle, method="most_probable")
        elif not np.isnan(vTh):
            # check units of thermal velocity
            vTh = vTh.to(u.m / u.s)
    elif np.isnan(vTh) and units == "unitless":
        # assuming unitless temperature is in Kelvins
        vTh = (
            parameters.thermal_speed(T * u.K, particle=particle, method="most_probable")
        ).si.value
    # Get thermal velocity squared
    vThSq = vTh ** 2
    # Get square of relative particle velocity
    vSq = (v - v_drift) ** 2
    # calculating distribution function
    coeff = (vThSq * np.pi) ** (-1 / 2)
    expTerm = np.exp(-vSq / vThSq)
    distFunc = coeff * expTerm
    if units == "units":
        return distFunc.to(u.s / u.m)
    elif units == "unitless":
        return distFunc

>>> Maxwellian_speed_1D(v=v, T=30000 * u.K, particle='e', v_drift=0 * u.m / u.s)
    

    """
    if units == "units":
        # unit checks and conversions
        # checking velocity units
        v = v.to(u.m / u.s)
        # Catching case where drift velocities have default values, they
        # need to be assigned units
        v_drift = _v_drift_units(v_drift)
        # convert temperature to Kelvins
        T = T.to(u.K, equivalencies=u.temperature_energy())
        if np.isnan(vTh):
            # get thermal velocity and thermal velocity squared
            vTh = parameters.thermal_speed(T, particle=particle, method="most_probable")
        elif not np.isnan(vTh):
            # check units of thermal velocity
            vTh = vTh.to(u.m / u.s)
    elif np.isnan(vTh) and units == "unitless":
        # assuming unitless temperature is in Kelvins
        vTh = (
            parameters.thermal_speed(T * u.K, particle=particle, method="most_probable")
        ).si.value
    # Get thermal velocity squared
    vThSq = vTh ** 2
    # Get square of relative particle velocity
    vSq = (v - v_drift) ** 2
    # calculating distribution function
    coeff = 2 * (vThSq * np.pi) ** (-1 / 2)
    expTerm = np.exp(-vSq / vThSq)
    distFunc = coeff * expTerm

"""
    if v_drift != 0:
        raise NotImplementedError("Non-zero drift speed is work in progress.")
    if units == "units":
        # unit checks and conversions
        # checking velocity units
        v = v.to(u.m / u.s)
        # Catching case where drift velocity has default value, and
        # needs to be assigned units
        v_drift = _v_drift_units(v_drift)
        # convert temperature to Kelvins
        T = T.to(u.K, equivalencies=u.temperature_energy())
        if np.isnan(vTh):
            # get thermal velocity and thermal velocity squared
            vTh = parameters.thermal_speed(T, particle=particle, method="most_probable")
        elif not np.isnan(vTh):
            # check units of thermal velocity
            vTh = vTh.to(u.m / u.s)
    elif np.isnan(vTh) and units == "unitless":
        # assuming unitless temperature is in Kelvins
        vTh = (
            parameters.thermal_speed(T * u.K, particle=particle, method="most_probable")
        ).si.value
    # getting square of thermal speed
    vThSq = vTh ** 2
    # get square of relative particle speed
    vSq = (v - v_drift) ** 2
    # calculating distribution function
    coeff1 = (np.pi * vThSq) ** (-1)
    coeff2 = 2 * np.pi * (v - v_drift)
    expTerm = np.exp(-vSq / vThSq)

-------
    >>> from astropy import units as u
    >>> from numpy import pi
    >>> from astropy.constants import c
    >>> T = 30 * 11600 * u.K
    >>> n = 1e18 * u.cm**-3
    >>> particle = 'Ne'
    >>> z_mean = 8 * u.dimensionless_unscaled
    >>> vTh = parameters.thermal_speed(T, particle, method="most_probable")
    >>> omega = 5.635e14 * 2 * pi * u.rad / u.s
    >>> kWave = omega / vTh
    >>> permittivity_1D_Maxwellian(omega, kWave, T, n, particle, z_mean)
    
    """
    # thermal velocity
    vTh = parameters.thermal_speed(T=T, particle=particle, method="most_probable")
    # plasma frequency
    wp = parameters.plasma_frequency(n=n, particle=particle, z_mean=z_mean)
    # scattering parameter alpha.
    # explicitly removing factor of sqrt(2) to be consistent with Froula
    alpha = np.sqrt(2) * (wp / (kWave * vTh)).to(u.dimensionless_unscaled)
    # The dimensionless phase velocity of the propagating EM wave.
    zeta = (omega / (kWave * vTh)).to(u.dimensionless_unscaled)
    chi = alpha ** 2 * (-1 / 2) * plasma_dispersion_func_deriv(zeta.value)
    return chi

if units == "units":
        # unit checks and conversions
        # checking velocity units
        vx = vx.to(u.m / u.s)
        vy = vy.to(u.m / u.s)
        vz = vz.to(u.m / u.s)
        # catching case where drift velocities have default values, they
        # need to be assigned units
        vx_drift = _v_drift_units(vx_drift)
        vy_drift = _v_drift_units(vy_drift)
        vz_drift = _v_drift_units(vz_drift)
        # convert temperature to Kelvins
        T = T.to(u.K, equivalencies=u.temperature_energy())
        if np.isnan(vTh):
            # get thermal velocity and thermal velocity squared
            vTh = parameters.thermal_speed(T, particle=particle, method="most_probable")
        elif not np.isnan(vTh):
            # check units of thermal velocity
            vTh = vTh.to(u.m / u.s)
    elif np.isnan(vTh) and units == "unitless":
        # assuming unitless temperature is in Kelvins
        vTh = parameters.thermal_speed(
            T * u.K, particle=particle, method="most_probable"
        ).si.value
    # accounting for thermal velocity in 3D
    vThSq = vTh ** 2
    # Get square of relative particle velocity
    vSq = (vx - vx_drift) ** 2 + (vy - vy_drift) ** 2 + (vz - vz_drift) ** 2
    # calculating distribution function
    coeff = (vThSq * np.pi) ** (-3 / 2)
    expTerm = np.exp(-vSq / vThSq)
    distFunc = coeff * expTerm

def _replaceNanVwithThermalV(V, T, m):
    """
    Get thermal velocity of system if no velocity is given, for a given mass.
    Handles vector checks for V, you must already know that T and m are okay.
    """
    if np.any(V == 0):
        raise utils.PhysicsError("You cannot have a collision for zero velocity!")
    # getting thermal velocity of system if no velocity is given
    if V is None:
        V = parameters.thermal_speed(T, mass=m)
    elif np.any(np.isnan(V)):
        if np.isscalar(V.value) and np.isscalar(T.value):
            V = parameters.thermal_speed(T, mass=m)
        elif np.isscalar(V.value):
            V = parameters.thermal_speed(T, mass=m)
        elif np.isscalar(T.value):
            V = V.copy()
            V[np.isnan(V)] = parameters.thermal_speed(T, mass=m)
        else:
            V = V.copy()
            V[np.isnan(V)] = parameters.thermal_speed(T[np.isnan(V)], mass=m)
    return V