How to use the pylops.signalprocessing.Convolve1D function in pylops

np.diag(G3[itheta] * np.ones(nt0))))
             for itheta in range(ntheta)]
        G = np.vstack(G).reshape(ntheta * nt0, 3 * nt0)

        # Create wavelet operator
        C = convmtx(wav, nt0)[:, len(wav) // 2:-len(wav) // 2 + 1]
        C = [C] * ntheta
        C = block_diag(*C)

        # Combine operators
        M = np.dot(C, np.dot(G, D))
        return MatrixMult(M, dims=spatdims)

    else:
        # Create wavelet operator
        Cop = Convolve1D(np.prod(np.array(dims)), h=wav,
                         offset=len(wav)//2, dir=0, dims=dims)

        # create AVO operator
        AVOop = AVOLinearModelling(theta, vsvp, spatdims=spatdims,
                                   linearization=linearization)

        # Create derivative operator
        dimsm = list(dims)
        dimsm[1] = AVOop.npars
        Dop = FirstDerivative(np.prod(np.array(dimsm)), dims=dimsm,
                              dir=0, sampling=1.)
        return Cop*AVOop*Dop

"""
    np.random.seed(10)
    # 1D
    if par['dir'] == 0:
        Cop = Convolve1D(par['nx'], h=h1, offset=par['offset'],
                         dtype='float32')
        assert dottest(Cop, par['nx'], par['nx'])

        x = np.zeros((par['nx']))
        x[par['nx']//2] = 1.
        xlsqr = lsqr(Cop, Cop * x, damp=1e-20, iter_lim=200, show=0)[0]
        assert_array_almost_equal(x, xlsqr, decimal=1)

    # 1D on 2D
    if par['dir'] < 2:
        Cop = Convolve1D(par['ny'] * par['nx'], h=h1, offset=par['offset'],
                        dims=(par['ny'], par['nx']), dir=par['dir'],
                        dtype='float32')
        assert dottest(Cop, par['ny'] * par['nx'], par['ny'] * par['nx'])

        x = np.zeros((par['ny'], par['nx']))
        x[int(par['ny']/2-3):int(par['ny']/2+3),
          int(par['nx']/2-3):int(par['nx']/2+3)] = 1.
        x = x.flatten()
        xlsqr = lsqr(Cop, Cop * x, damp=1e-20, iter_lim=200, show=0)[0]
        assert_array_almost_equal(x, xlsqr, decimal=1)

    # 1D on 3D
    Cop = Convolve1D(par['nz'] * par['ny'] * par['nx'], h=h1,
                     offset=par['offset'],
                     dims=(par['nz'], par['ny'], par['nx']), dir=par['dir'],
                     dtype='float32')

itrav = (trav / dt).astype('int32')
        travd = (trav / dt - itrav)
        if ndim == 2:
            itrav = itrav.reshape(nx, nz, ns * nr)
            travd = travd.reshape(nx, nz, ns * nr)
            dims = tuple(dims)
        else:
            itrav = itrav.reshape(ny*nx, nz, ns * nr)
            travd = travd.reshape(ny*nx, nz, ns * nr)
            dims = (dims[0]*dims[1], dims[2])

        # create operator
        sop = Spread(dims=dims, dimsd=(ns * nr, nt),
                     table=itrav, dtable=travd, engine='numba')

        cop = Convolve1D(ns * nr * nt, h=wav, offset=wavcenter,
                         dims=(ns * nr, nt),
                         dir=1)
        demop = cop * sop
    else:
        raise NotImplementedError('method must be analytic or eikonal')
    return demop

xinv = Cop / y

fig, ax = plt.subplots(1, 1, figsize=(10, 3))
ax.plot(t, x, 'k', lw=2, label=r'$x$')
ax.plot(t, y, 'r', lw=2, label=r'$y=Ax$')
ax.plot(t, xinv, '--g', lw=2, label=r'$x_{ext}$')
ax.set_title('Convolve 1d data', fontsize=14, fontweight='bold')
ax.legend()
ax.set_xlim(1.9, 2.1)

###############################################################################
# We show now that also a filter with mixed phase (i.e., not centered
# around zero) can be applied and inverted for using the
# :py:class:`pylops.signalprocessing.Convolve1D`
# operator.
Cop = pylops.signalprocessing.Convolve1D(nt, h=h, offset=hcenter - 3,
                                         dtype='float32')
y = Cop*x
y1 = Cop.H*x
xinv = Cop / y

fig, ax = plt.subplots(1, 1, figsize=(10, 3))
ax.plot(t, x, 'k', lw=2, label=r'$x$')
ax.plot(t, y, 'r', lw=2, label=r'$y=Ax$')
ax.plot(t, y1, 'b', lw=2, label=r'$y=A^Hx$')
ax.plot(t, xinv, '--g', lw=2, label=r'$x_{ext}$')
ax.set_title('Convolve 1d data with non-zero phase filter', fontsize=14,
             fontweight='bold')
ax.set_xlim(1.9, 2.1)
ax.legend()

###############################################################################

###############################################################################
# We will start by creating a zero signal of lenght :math:`nt` and we will
# place a unitary spike at its center. We also create our filter to be
# applied by means of :py:class:`pylops.signalprocessing.Convolve1D` operator.
# Following the seismic example mentioned above, the filter is a
# `Ricker wavelet `_
# with dominant frequency :math:`f_0 = 30 Hz`.
nt = 1001
dt = 0.004
t = np.arange(nt)*dt
x = np.zeros(nt)
x[int(nt/2)] = 1
h, th, hcenter = ricker(t[:101], f0=30)

Cop = pylops.signalprocessing.Convolve1D(nt, h=h, offset=hcenter,
                                         dtype='float32')
y = Cop*x

xinv = Cop / y

fig, ax = plt.subplots(1, 1, figsize=(10, 3))
ax.plot(t, x, 'k', lw=2, label=r'$x$')
ax.plot(t, y, 'r', lw=2, label=r'$y=Ax$')
ax.plot(t, xinv, '--g', lw=2, label=r'$x_{ext}$')
ax.set_title('Convolve 1d data', fontsize=14, fontweight='bold')
ax.legend()
ax.set_xlim(1.9, 2.1)

###############################################################################
# We show now that also a filter with mixed phase (i.e., not centered
# around zero) can be applied and inverted for using the

and the third direction:

    .. math::
        y[i,j,k] = 1/n_{smooth} \sum_{l=-(n_{smooth}-1)/2}^{(n_{smooth}-1)/2}
        x[i,j,l]

    Note that since the filter is symmetrical, the *Smoothing1D* operator is
    self-adjoint.

    """
    if isinstance(dims, int):
        dims = (dims,)
    if nsmooth % 2 == 0:
        nsmooth += 1

    return Convolve1D(np.prod(np.array(dims)), dims=dims, dir=dir,
                      h=np.ones(nsmooth)/float(nsmooth), offset=(nsmooth-1)/2,
                      dtype=dtype)

# :py:class:`pylops.optimization.sparsity.FISTA` solvers allows us
# to succesfully retrieve the input signal even in the presence of noise.
# :py:class:`pylops.optimization.sparsity.FISTA` shows faster convergence which
# is particularly useful for this problem.

nt = 61
dt = 0.004
t = np.arange(nt)*dt
x = np.zeros(nt)
x[10] = -.4
x[int(nt/2)] = 1
x[nt-20] = 0.5

h, th, hcenter = pylops.utils.wavelets.ricker(t[:101], f0=20)

Cop = pylops.signalprocessing.Convolve1D(nt, h=h, offset=hcenter,
                                         dtype='float32')
y = Cop*x
yn = y + np.random.normal(0, 0.1, y.shape)

# noise free
xls = Cop / y

xomp, nitero, costo = \
    pylops.optimization.sparsity.OMP(Cop, y, niter_outer=200, sigma=1e-8)

xista, niteri, costi = \
    pylops.optimization.sparsity.ISTA(Cop, y, niter=400, eps=5e-1,
                                      tol=1e-8, returninfo=True)

fig, ax = plt.subplots(1, 1, figsize=(8, 3))
ax.plot(t, x, 'k', lw=8, label=r'$x$')