home *** CD-ROM | disk | FTP | other *** search
- /* CMD: Surface Plot
- *
- * 2-D Function Surface Maker for Modeler
- * Originally by Arnie Cachelin © 1992 NewTek Inc., Sat Jul 11 1992
- *
- * Modified to use Modeler-based UI by Stuart Ferguson, 11/92
- * Further modified by Arnie, Sat May 8 16:04:25 1993
- */
-
- arg bas
- call addlib "LWModelerARexx.port", 0
- call addlib "rexxsupport.library", 0, -30, 0
- signal on error
- signal on syntax
-
- MATHLIB="rexxmathlib.library"
- IF POS(MATHLIB , SHOW('L')) = 0 THEN
- IF ~ADDLIB(MATHLIB , 0 , -30 , 0) THEN DO
- call notify(1,"!Can't find "MATHLIB)
- exit
- END
-
- sysnam = 'Plot 2D Function'
- filnam = 'ENV:plot2d.state'
- version = 'Plot 2D v1.1'
-
- /* Setup state. Read stored one, if any.
- */
- NSX = 20
- NSY = 20
- xmin = -10
- xmax = 10
- ymin = -10
- ymax = 10
- func = '3*sin(x*y)'
- ufunc=func
- fn = 3
- flip = 1
- tri = 1
- if bas="" then bas = 0
- bas = 0
- typ=2
- if (exists(filnam)) then do
- if (~open(state, filnam, 'R')) then break
- if (readln(state) ~= version) then break
- parse value readln(state) with nsx nsy xmin xmax ymin ymax fn tri typ .
- func = readln(state)
- call close state
- end
-
- FnList.1='RandomSheet'
- FnList.2='Wave'
- FnList.3='RadialWave'
- FnList.4='Gaussian'
- FnList.5='Interfere'
- FnList.6='Custom'
- FList= FnList.1 FnList.2 FnList.3 FnList.4 FnList.5 FnList.6
-
- call req_begin sysnam
-
- id_fnl = req_addcontrol("f(x,y)", 'CH',FList)
- id_lox = req_addcontrol("Low X", 'n', 1)
- id_hix = req_addcontrol("High X", 'n', 1)
- id_nsx = req_addcontrol("X Segments", 'n')
- id_loy = req_addcontrol("Low Y", 'n', 1)
- id_hiy = req_addcontrol("High Y", 'n', 1)
- id_nsy = req_addcontrol("Y Segments", 'n')
- id_fun = req_addcontrol("Custom Function", 's', 35)
- id_tri = req_addcontrol("Triangles", 'b')
- id_typ = req_addcontrol("Build: ","CH","Points Polys Curves")
-
- call req_setval id_lox, xmin, -10
- call req_setval id_loy, ymin, -10
- call req_setval id_hix, xmax, 10
- call req_setval id_hiy, ymax, 10
- call req_setval id_nsx, nsx, 20
- call req_setval id_nsy, nsy, 20
- call req_setval id_fun, func
- call req_setval id_fnl, fn,6
- call req_setval id_tri, tri,tri
- call req_setval id_typ, typ,typ
-
- if (~req_post()) then do
- call req_end
- exit
- end
-
- NSX = req_getval(id_nsx) % 1
- NSY = req_getval(id_nsx) % 1
- xmin = req_getval(id_lox)
- xmax = req_getval(id_hix)
- ymin = req_getval(id_loy)
- ymax = req_getval(id_hiy)
- fn = req_getval(id_fnl)
- typ = req_getval(id_typ)
- if fn=6 then func = req_getval(id_fun)
- else func = FnList.fn'(x,y)'
- tri = req_getval(id_tri)
-
- call req_end
-
- if (open(state, filnam, 'W')) then do
- call writeln state, version
- call writeln state, nsx nsy xmin xmax ymin ymax fn tri typ
- call writeln state, ufunc
- call close state
- end
-
- xrange = xmax - xmin
- yrange = ymax - ymin
- xmesh = xrange / NSX
- ymesh = yrange / NSY
- rscale = sqrt(xrange*yrange)
- tri_height = sqrt(3)/2
- ifunc = "z =" func
- say ifunc
- z=0
- zmax=z
- zmin=z
- call randu(time('s')) /* Seed random number generator */
-
- if typ=3 then call RectCurves
- else
- if tri then do
- if typ=1 then call TriPoints
- else call TriMesh
- end
- else do
- if typ=1 then call RectPoints
- else call RectMesh
- end
-
- /* if Bas totalpoints=totalpoints+MakeBase() */
-
- l1 = "Points created:" totalpoints
- l2 = "Polygons created:" poly
- l3 = "Z ranges between" zmin "and" zmax
- call notify 1, '!'sysnam, l1, l2, l3
-
- exit
-
- RectMesh:
- totalpoints = (NSX+1) * (NSY+1)
- totalpolys = NSX * NSY
- call add_begin
- call meter_begin totalpoints+2, sysnam, "Computing "totalpoints" points for "totalpolys" squares"
- do y=ymin to ymax by ymesh
- if y=ymax then TopCorner.4=totalpoints
- do x=xmin to xmax by xmesh
- interpret ifunc
- if z<zmin then zmin=z
- if z>zmax then zmax=z /* Just some silly stats for later */
- if (flip) then vec =x z y
- else vec = x y z
- call add_point(vec)
- call meter_step
- end
- if y=ymin then TopCorner.2=totalpoints
- if y=ymax then TopCorner.3=totalpoints
- end
-
- point=1
- poly=0
- call meter_begin totalpolys, sysnam, "Generating "totalpolys" Polygon Mesh"
- do y=ymin to ymax-ymesh by ymesh /* Don't make wrap-around polygon */
- do x=xmin to xmax by xmesh
- if x<xmax then do /* Again, Don't make wrap-around polygon! */
- if (flip) then
- call add_quad point point+NSX+1 point+NSX+2 point+1
- else
- call add_quad point point+1 point+NSX+2 point+NSX+1
- poly = poly + 1
- call meter_step
- end
- point = point + 1
- end
- end
- call meter_end
- call add_end
- return totalpoints
- /* */
-
- RectPoints:
- totalpoints = (NSX+1) * (NSY+1)
- poly=totalpoints
- call add_begin
- call meter_begin totalpoints+2, sysnam, "Computing "totalpoints" points for 1-point polygons"
- point=1
- do y=ymin to ymax by ymesh
- do x=xmin to xmax by xmesh
- interpret ifunc
- if z<zmin then zmin=z
- if z>zmax then zmax=z /* Just some silly stats for later */
- if (flip) then vec =x z y
- else vec = x y z
- call add_point(vec)
- call add_polygon(point)
- point=point+1
- call meter_step
- end
- end
- call add_end
- return totalpoints
- /* */
-
- RectCurves:
- totalpoints = (NSX+1) * (NSY+1)
- poly=NSX+1 + NSY+1
- call add_begin
- call meter_begin totalpoints+NSX+2, sysnam, "Computing "totalpoints" points for "poly" Curves"
- crv=""
- point=1
- do y=ymin to ymax by ymesh
- do x=xmin to xmax by xmesh
-
- interpret ifunc
-
- if z<zmin then zmin=z
- if z>zmax then zmax=z
-
- if (flip) then vec =x z y
- else vec = x y z
- call add_point(vec)
- crv=crv point
- point=point+1
- call meter_step
- end
- call Add_Curve(crv)
- crv=""
- end
- do p=1 to NSX+1
- do o=0 to NSY
- crv=crv p+o*(NSX+1)
- end
- call meter_step
- call Add_Curve(crv)
- crv=""
- end
- call add_end
- return totalpoints
- /* */
-
- TriMesh:
- totalpoints = (NSX+1) * (NSY+1)
- totalpolys = NSX * NSY * 2
- call add_begin
- call meter_begin totalpoints*2, sysnam, "Computing "totalpoints" points"
- offset=0
- rows=0
- totalpoints=0
-
- do y=ymin to ymax by ymesh* tri_height
- rows=rows+1
- columns=0
- if y=ymax then TopCorner.4=totalpoints
- do x=xmin+offset to xmax+offset by xmesh
-
- columns = columns + 1
- interpret ifunc
-
- if z<zmin then zmin=z
- if z>zmax then zmax=z /* Just some silly stats for later */
-
- if (flip) then vec =x z y
- else vec = x y z
- call add_point vec
- totalpoints = totalpoints + 1
- call meter_step
- end
- if y=ymin then TopCorner.2=totalpoints
- if y=ymax then TopCorner.3=totalpoints
- if offset=0 then offset=.5 * xmesh /* offset alternate lines */
- else offset=0
- end
- call meter_end
-
- point=1
- poly=0
- off=0
- call meter_begin totalpolys, sysnam, "Generating "totalpolys" Polygon Mesh"
- do row=0 to rows-2
- if off=0 then off=1 /* Boy this feels kludgey!!! */
- else off=0
- do col=1 to columns - 1
- if (flip) then do
- call add_quad col+row*columns col+(row*columns)+1 col+((row+1)*columns)+abs(off-1)
- call add_quad col+(row*columns)+off col+((row+1)*columns)+1 col+((row+1)*columns)
- poly=poly+2
- end
- else do
- call add_quad col+row*columns col+((row+1)*columns)+abs(off-1) col+(row*columns)+1
- call add_quad col+(row*columns)+off col+((row+1)*columns) col+((row+1)*columns)+1
- poly=poly+2
- end
- call meter_step
- end
- end
- call meter_end
- call add_end
- return totalpoints
- /* */
-
- TriPoints:
- totalpoints = (NSX+1) * (NSY+1)
- call add_begin
- call meter_begin totalpoints*2, sysnam, "Computing "totalpoints" points for 1-point polygons"
- offset=0
- rows=0
- totalpoints=0
- point=1
-
- do y=ymin to ymax by ymesh* tri_height
- rows=rows+1
- columns=0
- do x=xmin+offset to xmax+offset by xmesh
-
- columns = columns + 1
- interpret ifunc
-
- if z<zmin then zmin=z
- if z>zmax then zmax=z /* Just some silly stats for later */
-
- if (flip) then vec =x z y
- else vec = x y z
- call add_point vec
- call add_polygon point
- point=point+1
- call meter_step
- end
- if offset=0 then offset=.5 * xmesh /* offset alternate lines */
- else offset=0
- end
- call meter_end
- totalpoints = point - 1
- poly=totalpoints
- call add_end()
- return totalpoints
- /* */
-
-
- MakeBase:
- add_begin
- point=TotalPoints-1
- say point
- TopCorner.1=1
- BotCorner.1=totalpoints
- BotCorner.2=BotCorner.1 + NSX +1
- BotCorner.3=BotCorner.2 + NSY
- BotCorner.4=BotCorner.3 + NSX +1
- z=ZMin-(xrange+yrange)/4
- y=ymin
- BotCorner.1=point+1
- do x=xmin to xmax by xmesh
- call add_point x y z
- point=point+1
- pointList.point = x y z
- end
- BotCorner.2=point
- x=xmin
- do y=ymin to ymax by ymesh
- call add_point x y z
- point=point+1
- pointList.point = x y z
- end
- y=ymax
- BotCorner.4=point+1
- do x=xmin to xmax by xmesh
- call add_point x y z
- point=point+1
- pointList.point = x y z
- end
- BotCorner.3=point
- x=xmax
- do y=ymin to ymax by ymesh
- call add_point x y z
- point=point+1
- pointList.point = x y z
- end
-
- call surface("AreaBottom")
- call add_quad BotCorner.1 BotCorner.2 BotCorner.3 BotCorner.4
- poly=poly+1
- call surface("AreaBase")
-
- do i=0 to NSX-1
- call add_quad TopCorner.1+i TopCorner.1+i+1 BotCorner.1+i+1 BotCorner.1+i
- poly=poly+1
- end
-
- do i=0 to NSY-1
- call add_quad TopCorner.4+i TopCorner.4+i+1 BotCorner.4+i+1 BotCorner.4+i
- poly=poly+1
- end
-
- do i=0 to NSX-1
- a=((TopCorner.1)+i*(NSX+1))
- b=((TopCorner.1)+(i+1)*(NSX+1))
- c=((BotCorner.2)+i+1)
- d=((BotCorner.2)+i+2)
- call add_quad a b d c
- poly=poly+1
- end
-
- do i=0 to NSY-1
- a=((TopCorner.2)+i*(NSY+1))
- b=((TopCorner.2)+(i+1)*(NSY+1))
- c=((BotCorner.3)+i+1)
- d=((BotCorner.3)+i+2)
- call add_quad a b d c
- poly=poly+1
- end
-
- add_end
- return point
- /* */
-
- Radius: PROCEDURE
- arg xf, yf
- return sqrt(xf*xf+yf*yf)
-
- Sinc: PROCEDURE EXPOSE rscale
- arg xf, yf /* Classic the spherical Bessel f'n j0 */
- r=Radius(xf,yf)
- if (r = 0)
- zf = rscale
- else
- zf = rscale * sin(r) / r
- return zf
-
- Wave: PROCEDURE EXPOSE rscale
- arg xf, yf /* Simple wavy sheet */
- zf=rscale*sin(2*3.141592*xf/rscale)/3
- return(zf)
-
- RadialWave: PROCEDURE EXPOSE rscale
- arg xf, yf /* rings of waves sheet */
- xc=0; yc=0 /* Center coord.s */
- r=Radius(xf-xc,yf-yc)
- zf=rscale*sin(10*3.141592*r/rscale)/5
- return(zf)
-
- Interfere: PROCEDURE EXPOSE rscale
- arg xf, yf /* Interference of several (3) sources */
- r=Radius(xf,yf) /* (0,0) */
- d1=Radius(xf-5*rscale/10,yf+2*rscale/10) /* (5,-2) */
- d2=Radius(xf+4*rscale/10,yf-3*rscale/10) /* (-4,3) */
- zf=(rscale/7)*(sin(r)+sin(d1*1.5)+1.4*cos(d2))
- return(zf)
-
- RandomSheet: PROCEDURE EXPOSE rscale
- arg xf, yf /* Random altitudes */
- amp=(rscale/3)
- zf=randu()*amp-amp/2
- return(zf)
-
- Polynomial: PROCEDURE EXPOSE rscale
- arg xf, yf
- zf = 2*xf*xf + 3*xf - 2*yf*yf + 8*yf + 2
- zf = zf / 20
- return(zf)
-
- Gaussian: PROCEDURE EXPOSE rscale
- arg xf, yf
- xc=0; yc=0 /* Center coord.s */
- amp = rscale/3
- r = Radius(xf-xc,yf-yc)
- zf = amp * exp(-r*r/20)
- return(zf)
-
- syntax:
- error:
- call end_all
- t=Notify(1,'!Rexx Script Error','@'ErrorText(rc),'Line 'SIGL)
- exit