; file vectorsurf.pro
;--------------------------------------------------------------------
;+
;  FILE:
;       vectorsurf.pro
;  CALLING SEQUENCE: vectorsurf              ; demo mode
;                    vectorsurf, x, y, z     ; 1-d vectors of same length
;                    vectorsurf, x, y, z, bg_color=100   ; chose background color
;                    vectorsurf, x, y, z, DRAWXSIZE = 600, DRAWYSIZE = 400
;  PURPOSE:
;       Shows a fancy way to plot irregulary spaced data from x,y & z vectors.
;
;  CATEGORY:
;       Surface Plotting
;  EXAMPLE:
;     IDL> x=[4,6,1,7,1,8,5,2,4]
;     IDL> y=[3,3,1,6,4,5,6,8,1]
;     IDL> z=sqrt((x-3)^2+(y-2)^2)+randomn(seed)+2
;     IDL> vectorsurf, x, y, z
;  MODIFICATION HISTORY:
;       05-Sep-97 Written by Bill Davis
;                 from IDL's d_mathstat.pro in the IDL 5.0 release
;                 Written by:  DC, RSI,  1995
;-
;--------------------------------------------------------------------

PRO vectorsurf, x, y, z, DRAWXSIZE = drawXsize, DRAWYSIZE = drawYsize, $
                BG_COLOR = bg_color

IF ( (!D.NAME NE 'X') AND (!D.NAME NE 'MAC') ) THEN BEGIN
    PRINT,'*** vectorsurf only works on X windows ***'
    RETURN
ENDIF

IF N_ELEMENTS(drawXsize) LE 0 THEN drawXsize = 400
IF N_ELEMENTS(drawYsize) LE 0 THEN drawYsize = 400

window,11, xsize = drawXsize, ysize = drawYsize

IF N_PARAMS(0) LE 0 THEN BEGIN
    PRINT, '> No parameters entered, will plot some demo data'
    PRINT, ' '
    PRINT, '  Normal call:'
    PRINT, '     vectorsurf, x, y, z'
    PRINT, ' '
    PRINT, '     Where x, y & z are 1-dimensional vectors of the same length'
    PRINT, ' '

    points=20

    x = RANDOMN(s, points)
    y = RANDOMN(s, points)
    z = SIN(x - COS(y)^2) - COS(SIN(x^2) + SIN(y)) + 4.0
ENDIF  ELSE points = N_ELEMENTS(z)

zz = MIN_CURVE_SURF(z, x, y, Nx=points, Ny=points)
min_x = MIN(x, MAX=max_x)
min_y = MIN(y, MAX=max_y)
xx = min_x + ((max_x - min_x) * FINDGEN(points) / FLOAT(points-1))
yy = min_y + ((max_y - min_y) * FINDGEN(points) / FLOAT(points-1))

save_name = !D.Name
SET_PLOT, 'Z'
DEVICE, Set_Resolution=[drawXSize, drawYSize]

IF N_ELEMENTS(bg_color) GT 0 THEN BEGIN
        SHADE_SURF, zz, xx, yy, XSTYLE=1, YSTYLE=1, ZSTYLE=1, $
            TICKLEN=(-0.02), AX=65, AZ=30, BACKGROUND=bg_color, SKIRT=0.0
            ;;;  Include: "SHADES = BYTSCL(zz)" for shading proportional to z value
ENDIF ELSE BEGIN
        SHADE_SURF, zz, xx, yy, XSTYLE=1, YSTYLE=1, ZSTYLE=1, $
            TICKLEN=(-0.02), AX=65, AZ=30, BACKGROUND=255, SKIRT=0.0
            ;;;  Include: "SHADES = BYTSCL(zz)" for shading proportional to z value
ENDELSE

img = TVRD(0, 0, drawXSize, drawYSize)
ERASE
TV, img

SURFACE, zz, xx, yy, /NODATA, /NOERASE, $
XSTYLE=1, YSTYLE=1, ZSTYLE=1, $
    TICKLEN=(-0.02), COLOR=3, AX=65, AZ=30, SKIRT=0.0, /SAVE

for i=0, (N_ELEMENTS(x)-1L) do begin
    PLOTS, x(i), y(i), 0.0, PSYM=1, /T3D, /Data, COLOR=3, $
        SYMSIZE=3.0
    PLOTS, [x(i), x(i)], [y(i), y(i)], [0.0, z(i)], $
        /T3D, /DATA, COLOR=2
endfor

SURFACE, zz, xx, yy, /SAVE, XSTYLE=5, YSTYLE=5, ZSTYLE=5, $
    TICKLEN=(-0.02), /NOERASE, COLOR=15, AX=65, AZ=30, SKIRT=0.0

EMPTY
img = TVRD(0, 0, drawXSize, drawYSize)
DEVICE, /CLOSE
SET_PLOT, save_name
TV, img
EMPTY

delt_x = 0.01 / !X.S(1)
delt_y = 0.01 / !Y.S(1)
for i=0, (N_ELEMENTS(x)-1L) do begin
   px = [x(i)-delt_x, x(i)+delt_x, x(i)+delt_x, x(i)-delt_x, x(i)-delt_x]
   py = [y(i)-delt_y, y(i)-delt_y, y(i)+delt_y, y(i)+delt_y, y(i)-delt_y]
   ix = FLOAT(points-1) * (px - min_x) / (max_x - min_x)
   iy = FLOAT(points-1) * (py - min_y) / (max_y - min_y)
   pz = INTERPOLATE(zz, ix, iy)
   PLOTS, px, py, pz, /T3D, /DATA, COLOR=0
endfor

STOP

EMPTY
RETURN
end