<SYSTEM_TASK:>
Solve the following problem using Python, implementing the functions described below, one line at a time
<END_TASK>
<USER_TASK:>
Description:
 def gpr(src, rec, depth, res, freqtime, cf, gain=None, ab=11, aniso=None, epermH=None, epermV=None, mpermH=None, mpermV=None, xdirect=False, ht='quad', htarg=None, ft='fft', ftarg=None, opt=None, loop=None, verb=2):
 r"""Return the Ground-Penetrating Radar signal. THIS FUNCTION IS EXPERIMENTAL, USE WITH CAUTION. It is rather an example how you can calculate GPR responses; however, DO NOT RELY ON IT! It works only well with QUAD or QWE (``quad``, ``qwe``) for the Hankel transform, and with FFT (``fft``) for the Fourier transform. It calls internally ``dipole`` for the frequency-domain calculation. It subsequently convolves the response with a Ricker wavelet with central frequency ``cf``. If signal!=None, it carries out the Fourier transform and applies a gain to the response. For input parameters see the function ``dipole``, except for: Parameters cf : float Centre frequency of GPR-signal, in Hz. Sensible values are between 10 MHz and 3000 MHz. gain : float Power of gain function. If None, no gain is applied. Only used if signal!=None. Returns ------- EM : ndarray GPR response """

    if verb > 2:
        print("   GPR             :  EXPERIMENTAL, USE WITH CAUTION")
        print("     > centre freq :  " + str(cf))
        print("     > gain        :  " + str(gain))

    # === 1.  CHECK TIME ============

    # Check times and Fourier Transform arguments, get required frequencies
    time, freq, ft, ftarg = check_time(freqtime, 0, ft, ftarg, verb)

    # === 2. CALL DIPOLE ============

    EM = dipole(src, rec, depth, res, freq, None, ab, aniso, epermH, epermV,
                mpermH, mpermV, xdirect, ht, htarg, ft, ftarg, opt, loop, verb)

    # === 3. GPR STUFF

    # Get required parameters
    src, nsrc = check_dipole(src, 'src', 0)
    rec, nrec = check_dipole(rec, 'rec', 0)
    off, _ = get_off_ang(src, rec, nsrc, nrec, 0)

    # Reshape output from dipole
    EM = EM.reshape((-1, nrec*nsrc), order='F')

    # Multiply with ricker wavelet
    cfc = -(np.r_[0, freq[:-1]]/cf)**2
    fwave = cfc*np.exp(cfc)
    EM *= fwave[:, None]

    # Do f->t transform
    EM, conv = tem(EM, off, freq, time, 0, ft, ftarg)

    # In case of QWE/QUAD, print Warning if not converged
    conv_warning(conv, ftarg, 'Fourier', verb)

    # Apply gain; make pure real
    EM *= (1 + np.abs((time*10**9)**gain))[:, None]
    EM = EM.real

    # Reshape for number of sources
    EM = np.squeeze(EM.reshape((-1, nrec, nsrc), order='F'))

    return EM