<SYSTEM_TASK:>
Solve the following problem using Python, implementing the functions described below, one line at a time
<END_TASK>
<USER_TASK:>
Description:
 def burst_range(psd, snr=8, energy=1e-2, fmin=100, fmax=500):
 """Calculate the integrated GRB-like GW burst range from a strain PSD Parameters psd : `~gwpy.frequencyseries.FrequencySeries` the instrumental power-spectral-density data snr : `float`, optional the signal-to-noise ratio for which to calculate range, default: ``8`` energy : `float`, optional the relative energy output of the GW burst, defaults to ``1e-2`` for a GRB-like burst fmin : `float`, optional the lower frequency cutoff of the burst range integral, default: ``100 Hz`` fmax : `float`, optional the upper frequency cutoff of the burst range integral, default: ``500 Hz`` Returns ------- range : `~astropy.units.Quantity` the GRB-like-burst sensitive range [Mpc (default)] Examples -------- Grab some data for LIGO-Livingston around GW150914 and generate a PSD Now we can calculate the :func:`burst_range`: 42.5055584195 Mpc """

    freqs = psd.frequencies.value
    # restrict integral
    if not fmin:
        fmin = psd.f0
    if not fmax:
        fmax = psd.span[1]
    condition = (freqs >= fmin) & (freqs < fmax)
    # calculate integrand and integrate
    integrand = burst_range_spectrum(
        psd[condition], snr=snr, energy=energy) ** 3
    result = integrate.trapz(integrand.value, freqs[condition])
    # normalize and return
    r = units.Quantity(result / (fmax - fmin), unit=integrand.unit) ** (1/3.)
    return r.to('Mpc')