#
# smoothing operator demo
#
try:    from rsf.cluster import *
except: from rsf.proj    import *
import wplot
from smoLOP import *

from DPTEST import *
from LICGsolver import *
from LICGutils  import *

# ------------------------------------------------------------
do3D = 'y'
doINV = 'y'
# ------------------------------------------------------------

par = dict(
    nx=25,  ox=0, dx=1.0,  lx='x', ux='',
    ny=11,  oy=0, dy=1.0,  ly='y', uy='',
    nz=15,  oz=0, dz=1.0,  lz='z', uz=''
    )
wplot.param(par)

# ------------------------------------------------------------
Flow('s2D',None,
    '''
    spike nsp=1 mag=1 k1=8 k2=13
    n1=%(nz)d o1=%(oz)g d1=%(dz)g label1="" unit1=""
    n2=%(nx)d o2=%(ox)g d2=%(dx)g label2="" unit2=""
    '''%par)
Result('s2D',wplot.igrey2d('',par))

if(do3D=='y'):
    Flow('s3D',None,
    '''
    spike nsp=1 mag=1 k1=8 k2=13 k3=6
    n1=%(nz)d o1=%(oz)g d1=%(dz)g label1="" unit1=""
    n2=%(nx)d o2=%(ox)g d2=%(dx)g label2="" unit2=""
    n3=%(ny)d o3=%(oy)g d3=%(dy)g label3="" unit3=""
    '''%par)
    Result('s3D','byte gainpanel=a pclip=100 |'+wplot.igrey3d('',par))

# ------------------------------------------------------------
# initialize operators
# ------------------------------------------------------------
S2D = smo2d(3,3)
if(do3D=='y'):
    S3D = smo3d(3,3,3)

# ------------------------------------------------------------
# FORWARD/ADJOINT operators
# ------------------------------------------------------------
S2D.FORW('s2D',   's2Dfor')
S2D.ADJT('s2Dadj','s2Dfor')
Result('s2Dfor',wplot.igrey2d('',par))
Result('s2Dadj',wplot.igrey2d('',par))

if(do3D=='y'):
    S3D.FORW('s3D',   's3Dfor')
    S3D.ADJT('s3Dadj','s3Dfor')
    Result('s3Dfor','byte gainpanel=a pclip=100 |'+wplot.igrey3d('',par))
    Result('s3Dadj','byte gainpanel=a pclip=100 |'+wplot.igrey3d('',par))
    
# ------------------------------------------------------------
# DOT PRODUCT TEST
# ------------------------------------------------------------
Flow('S2Dd','s2Dfor','window squeeze=n')
Flow('S2Dm','s2Dadj','window squueze=n')
DP2D = DPTEST(S2D,['S2Dm'],['S2Dd'])
DP2D.RUN()

if(do3D=='y'):
    Flow('S3Dd','s3Dfor','window squeeze=n')
    Flow('S3Dm','s3Dadj','window squueze=n')
    DP3D = DPTEST(S3D,['S3Dm'],['S3Dd'])
    DP3D.RUN()

# ------------------------------------------------------------

if(doINV == 'y'):
    NITER = 21
    
    Flow('z2D','s2D','math output=0')

    # OLD
    OSMO = LICGOLD('OLD',par,
                       Lop=S2D,
                       d='s2Dfor',
                       xo='z2D',
                       xf='s2Dold',
                       niter=NITER)
    OSMO.iterate()
    plotOFUN('OLD',NITER,par)
    Result('s2Dold',wplot.igrey2d('',par))

    # NEW
    NSMO = LICG('NEW',par,
                    Lop = S2D,
                    d   = ['s2Dfor'],
                    xo  = ['z2D'],
                    xf  = ['s2Dnew'], niter=NITER)
    NSMO.iterate()
    plotOFUN('NEW',NITER,par)
    Result('s2Dnew',wplot.igrey2d('',par))
    
End()
