home *** CD-ROM | disk | FTP | other *** search
/ ProfitPress Mega CDROM2 …eeware (MSDOS)(1992)(Eng) / ProfitPress-MegaCDROM2.B6I / GRAPHICS / PLOT / 3DF.ZIP / 3DF.PAS < prev    next >
Encoding:
Pascal/Delphi Source File  |  1991-04-20  |  13.7 KB  |  472 lines
  1. { 3dF --- Simple viewer for real functions of form F(x,y).  }
  2. { Parse an expression to form a syntax tree.  Walk the tree }
  3. { to evaluate the expression on a grid over the x-y plane.  }
  4. { Draw orthographic projection with hidden lines removed,   }
  5. { view angle, etc. adjustable by user with cursor controls. }
  6. { Disclaimer:  This code is a mess.  It was written in a    }
  7. { very short time to prove that expensive hardware was not  }
  8. { required for a specific purpose.  Let the user beware.    }
  9.  
  10. { Compile with TP 5.0. }
  11.  
  12. Program p3dF;
  13. Uses Crt, Parse;
  14.  
  15. { Assembly language assist.  According to tprof, most
  16.   opportunity for getting faster is in DrawLine,
  17.   especially Plot (on my 4.77Mhz 8088 w/o NCP). }
  18. procedure SetVideo(mode: Byte); external;
  19. procedure ClearScreen; external;
  20. procedure SetPallette(pal: Byte); external;
  21. procedure SetBackGround(color: Byte); external;
  22. procedure MovePen(x, y: Word); external;
  23. procedure DrawLine(x, y: Word); external;
  24. procedure InitCloud; external;
  25. procedure UpdateCloud; external;
  26. procedure ClearImage; external;
  27. procedure ShowImage; external;
  28. {$ifdef profile}
  29. procedure initx; external;
  30. procedure inity; external;
  31. procedure initXinc; external;
  32. procedure initYinc; external;
  33. procedure plot; external;
  34. {$endif}
  35. {$L 3df.obj}
  36.  
  37. const
  38.   nGrids = 15;      { make me 10 for faster displays of fewer grids }
  39.   DevP = 192;       { vertical device size in pixels }
  40.   StepsInQuad = 18;
  41.   IntScale = 32;
  42.   ScaleBits = 5;
  43.   ScaleBitsM1 = 4;  { scale bits minus 1 (accounts for 2/1 aspect of pixels) }
  44.   pi = 3.1415926535897932385;  { Save fun call time. }
  45. var
  46.   eData: array[-nGrids..nGrids, -nGrids..nGrids] of Real;
  47.   eMin, eMax: Real;
  48.   yData: array[0..3, -nGrids..nGrids, -nGrids..nGrids] of Integer;
  49.   scrImg: array[0..199, 0..79] of Byte;
  50.   baseXTbl,
  51.   baseYTbl,
  52.   deltaXiTbl,
  53.   deltaYiTbl,
  54.   deltaXjTbl,                                    { how much to inc x for step in j }
  55.   deltaYjTbl: array[0..StepsInQuad] of Integer;  { how much to inc y for step in j }
  56.  
  57. procedure FillGrid(size, exag: Real; phi: Integer);
  58. var
  59.   e, de, i, j: Integer;
  60.   sPhi, cPhi,
  61.   imageHt, imageWd,
  62.   unitsToDev, eScale,
  63.   theta, dTheta, sTheta, cTheta, sizeP, scaleFac: Real;
  64. begin
  65.   sPhi := Sin(phi*pi/180);
  66.   cPhi := Cos(phi*pi/180);
  67.   { 1.5 is really sqrt(2)*fudge to _ensure_ images lie in screen. }
  68.   imageHt := sPhi*size*1.5+cPhi*(eMax-eMin)*exag;
  69.   imageWd := size*1.5;
  70.   if imageHt > imageWd/1.6 then
  71.     unitsToDev := DevP/imageHt
  72.   else
  73.     unitsToDev := DevP*1.6/imageWd;
  74.   sizeP := size*unitsToDev;
  75.   eScale := exag*unitsToDev*cPhi*IntScale;
  76.   for i := -nGrids to nGrids do
  77.     for j := -nGrids to nGrids do begin
  78.       e := Round((eData[i,j]-eMin)*eScale);
  79.       yData[0, i, j] := e;
  80.       yData[1, j,-i] := e;
  81.       yData[2, -i,-j] := e;
  82.       yData[3, -j, i] := e;
  83.     end;
  84.   theta := 0;
  85.   dTheta := (pi/2)/StepsInQuad;
  86.   scaleFac := sizeP/(nGrids+nGrids)*IntScale;
  87.   for i := 0 to stepsInQuad-1 do begin
  88.     sTheta := scaleFac*Sin(theta);
  89.     cTheta := scaleFac*Cos(theta);
  90.     deltaXiTbl[i] := Round(sTheta);
  91.     deltaYiTbl[i] := Round(cTheta*sPhi);
  92.     deltaXjTbl[i] := Round(cTheta);
  93.     deltaYjTbl[i] := Round(sTheta*sPhi);
  94.     baseXtbl[i] := Round((sizeP*0.707*Cos(theta+pi/4)+160)*IntScale);
  95.     baseYtbl[i] := Round((DevP-sizeP*0.707*(1-Sin(theta+pi/4))*sPhi)*IntScale);
  96.     theta := theta+dTheta;
  97.   end;
  98. end;
  99.  
  100. procedure DrawIt(cx, cy, size: Real);
  101. const
  102.   initStep = 4;
  103.   initDstep = 1;
  104.   initQuad = 1;
  105.   initPhi = 20;
  106.   initExag = 1;
  107. var
  108.   exag: Real;
  109.   i, j, phi,
  110.   step, quad, dStep, x00, y00, x0, y0, x, y, tx, ty,
  111.   px, py0, py, lastPx, lastPy0, dxi, dyi, dxj, dyj: Integer;
  112.   lastCutX, lastCutY: array[-nGrids..nGrids] of Integer;
  113.   ch: Char;
  114.  
  115.   procedure WritePhi;
  116.   begin
  117.     GoToXY(19, 25);
  118.     Write(phi:2);
  119.   end;
  120.  
  121.   procedure WriteExag;
  122.   begin
  123.     GoToXY(36, 25);
  124.     Write(exag:3:1);
  125.   end;
  126.  
  127. begin
  128.   DirectVideo := False;
  129.   dStep := initDstep;
  130.   step := initStep;
  131.   quad := initQuad;
  132.   phi := initPhi;
  133.   exag := initExag;
  134.   FillGrid(size, exag, phi);
  135.   SetVideo(6);  { 640 x 200 }
  136.   ClearScreen;
  137.   GoToXY(1, 25);
  138.   Write(
  139.     'eXit '#27#26'rot '#24#25'elev(__) PgUp/Dn exag(___) ',
  140.     'x ',cx:0:2, ' y ', cy:0:2, ' sz ', size:0:2,
  141.     ' f ', eMin:0:2, '/', eMax:0:2);
  142.   WritePhi;
  143.   WriteExag;
  144.   repeat
  145.     ClearImage;
  146.     InitCloud;
  147.     x0 := baseXtbl[step];
  148.     x00 := x0;
  149.     y0 := baseYtbl[step];
  150.     y00 := y0;
  151.     dxi := deltaXiTbl[step];
  152.     dyi := deltaYiTbl[step];
  153.     dxj := deltaXjTbl[step];
  154.     dyj := deltaYjTbl[step];
  155.     { draw first cut }
  156.     y := y0;
  157.     x := x0;
  158.     px := x shr ScaleBitsM1;
  159.     py0 := y shr ScaleBits;
  160.     py := (y-yData[quad, -nGrids,-nGrids]) shr ScaleBits;
  161.     MovePen(px, py0);
  162.     DrawLine(px, py);
  163.     lastPx := px;
  164.     lastPy0 := py0;
  165.     lastCutX[-nGrids] := px;
  166.     lastCutY[-nGrids] := py;
  167.     for j := -nGrids+1 to nGrids do begin
  168.       x := x - dxj;
  169.       y := y - dyj;
  170.       px := x shr ScaleBitsM1;
  171.       py0 := y shr ScaleBits;
  172.       py := (y-yData[quad, -nGrids, j]) shr ScaleBits;
  173.       DrawLine(px, py);
  174.       DrawLine(px, py0);
  175.       MovePen(lastPx, lastPy0);
  176.       DrawLine(px, py0);
  177.       if j = 0 then
  178.         if Odd(quad) then begin
  179.           tx := x-dxi*3;
  180.           ty := y+dyi*3;
  181.           DrawLine(tx shr ScaleBitsM1, ty shr ScaleBits);
  182.           if quad = 1 then begin
  183.             DrawLine((tx-dxi+dxj) shr ScaleBitsM1,
  184.              (ty+dyi+dyj) shr ScaleBits);
  185.             MovePen(tx shr ScaleBitsM1, ty shr ScaleBits);
  186.             DrawLine((tx-dxi-dxj) shr ScaleBitsM1,
  187.              (ty+dyi-dyj) shr ScaleBits);
  188.           end
  189.           else begin
  190.             DrawLine((tx+dxi+dxj) shr ScaleBitsM1,
  191.              (ty-dyi+dyj) shr ScaleBits);
  192.             MovePen(tx shr ScaleBitsM1, ty shr ScaleBits);
  193.             DrawLine((tx+dxi-dxj) shr ScaleBitsM1,
  194.              (ty-dyi-dyj) shr ScaleBits);
  195.           end;
  196.         end;
  197.       MovePen(px, py);
  198.       lastPx := px;
  199.       lastPy0 := py0;
  200.       lastCutX[j] := px;
  201.       lastCutY[j] := py;
  202.     end;
  203.     UpdateCloud;
  204.     lastPx := x00 shr ScaleBitsM1;
  205.     lastPy0 := y00 shr ScaleBits;
  206.     for i := -nGrids+1 to nGrids do begin
  207.       x0 := x0+dxi;
  208.       y0 := y0-dyi;
  209.       y := y0;
  210.       x := x0;
  211.       px := x shr ScaleBitsM1;
  212.       py0 := y shr ScaleBits;
  213.       py := (y-yData[quad, i,-nGrids]) shr ScaleBits;
  214.       MovePen(lastPx, lastPy0);
  215.       DrawLine(px, py0);
  216.       lastPx := px;
  217.       lastPy0 := py0;
  218.       if (i = 0) and (quad and 1 = 0) then begin
  219.         tx := x+dxj*3;
  220.         ty := y+dyj*3;
  221.         MovePen(tx shr ScaleBitsM1, ty shr ScaleBits);
  222.         if quad = 0 then begin
  223.           DrawLine((tx+dxj-dxi) shr ScaleBitsM1,
  224.             (ty+dyj+dyi) shr ScaleBits);
  225.           MovePen(tx shr ScaleBitsM1, ty shr ScaleBits);
  226.           DrawLine((tx+dxj+dxi) shr ScaleBitsM1,
  227.             (ty+dyj-dyi) shr ScaleBits);
  228.         end
  229.         else begin
  230.           DrawLine((tx-dxj-dxi) shr ScaleBitsM1,
  231.             (ty-dyj+dyi) shr ScaleBits);
  232.           MovePen(tx shr ScaleBitsM1, ty shr ScaleBits);
  233.           DrawLine((tx-dxj+dxi) shr ScaleBitsM1,
  234.             (ty-dyj-dyi) shr ScaleBits);
  235.         end;
  236.         MovePen(tx shr ScaleBitsM1, ty shr ScaleBits);
  237.         DrawLine(px, py0);
  238.       end;
  239.       DrawLine(px, py);
  240.       DrawLine(lastCutX[-nGrids], lastCutY[-nGrids]);
  241.       UpdateCloud;
  242.       MovePen(px, py);
  243.       lastCutX[-nGrids] := px;
  244.       lastCutY[-nGrids] := py;
  245.       for j := -nGrids+1 to nGrids do begin
  246.         x := x - dxj;
  247.         y := y - dyj;
  248.         px := x shr ScaleBitsM1;
  249.         py := (y-yData[quad, i, j]) shr ScaleBits;
  250.         DrawLine(px, py);
  251.         DrawLine(lastCutX[j], lastCutY[j]);
  252.         MovePen(px, py);
  253.         lastCutX[j] := px;
  254.         lastCutY[j] := py;
  255.         UpdateCloud;
  256.       end;
  257.     end;
  258.     if KeyPressed then begin
  259.       ch := ReadKey;
  260.       if ch = #0 then begin
  261.         ch := ReadKey;
  262.         case ch of
  263.         #77: if dStep > -3 then         { right arrow }
  264.                dStep := dStep-1;
  265.         #75: if dStep < 3 then          { left arrow  }
  266.                dStep := dStep+1;
  267.         #73: if exag < 9.5 then begin  { pg up }
  268.                if exag < 1 then
  269.                  exag := exag+0.1
  270.                else
  271.                  exag := exag+0.5;
  272.                FillGrid(size, exag, phi);
  273.                WriteExag;
  274.              end;
  275.         #81: if exag > 0.1 then begin     { pg dn }
  276.                if exag > 1 then
  277.                  exag := exag-0.5
  278.                else
  279.                  exag := exag-0.1;
  280.                FillGrid(size, exag, phi);
  281.                WriteExag;
  282.              end;
  283.         #72: if phi < 60 then begin     { up arrow }
  284.                phi := phi + 5;
  285.                FillGrid(size, exag, phi);
  286.                WritePhi;
  287.              end;
  288.         #80: if phi > 5 then begin      { dn arrow }
  289.                phi := phi - 5;
  290.                FillGrid(size, exag, phi);
  291.                WritePhi;
  292.              end;
  293.         #71: dStep := initDstep;        { home }
  294.         #79: begin                      { end }
  295.                phi := initPhi;
  296.                exag := initExag;
  297.                FillGrid(size, exag, phi);
  298.                WriteExag;
  299.                WritePhi;
  300.              end;
  301.         end;
  302.       end;
  303.     end
  304.     else ch := ' ';
  305.     while KeyPressed do ch := ReadKey; { purge keyboard buffer }
  306.     step := step+dStep;
  307.     if step >= stepsInQuad then begin
  308.       step := step mod stepsInQuad;
  309.       quad := quad+1;
  310.       if quad > 3 then
  311.         quad := 0;
  312.     end;
  313.     if step < 0 then begin
  314.       step := stepsInQuad+step mod stepsInQuad;
  315.       quad := quad-1;
  316.       if quad < 0 then
  317.         quad := 3;
  318.     end;
  319.     ShowImage;
  320.   until ch in [#27, 'X', 'x'];
  321.   SetVideo(CO80);
  322. end;
  323.  
  324. procedure ComputeFun(var expr: String; cx, cy, size: Real);
  325. var
  326.   tr: Tree;
  327.   i, j, n: Integer;
  328.   x, y, halfSize, e: Real;
  329. begin
  330.   tr := ParseExpr(expr);
  331.   halfSize := size/2;
  332.   eMax := -1e30; eMin := 1e30;
  333.   n := (nGrids + nGrids + 1);
  334.   n := n*n;
  335.   Write('Computing F (');
  336.   for j := -nGrids to nGrids do begin
  337.     Write(n:4, ' left).'#8#8#8#8#8#8#8#8#8#8#8);
  338.     x := cX + j/nGrids*halfSize;
  339.     for i := -nGrids to nGrids do begin
  340.       y := cY + i/nGrids*halfSize;
  341.       e := Eval(tr, x, y);
  342.       if e < eMin then eMin := e;
  343.       if e > eMax then eMax := e;
  344.       eData[i,j] := e;
  345.     end;
  346.     Dec(n, nGrids+nGrids+1);
  347.   end;
  348.   WriteLn(n:4);
  349. end;
  350.  
  351. procedure Usage(msg: String);
  352. const
  353.   nLines= 4;
  354.   verbiage: array[1..nLines] of string = (
  355.   'Usage: 3df [[options] expression]'#13#10+
  356.   '  options:'#13#10+
  357.   '    /x<real>    -- x of domain center (default 0)'#13#10+
  358.   '    /y<real>    -- y of domain center (default 0)'#13#10+
  359.   '    /s<posReal> -- x and y size of domain (default 1)'#13#10+
  360.   '    /?          -- This help.'#13#10,
  361.   '  expression: Optionally "quoted" function of x and y.'#13#10+
  362.   '    Operators by precedence:'#13#10+
  363.   '      a^b - a to b power, b >= 0 or int < 0    right associative'#13#10+
  364.   '      -a  - negate a'#13#10+
  365.   '      a*b - a times b    a/b - a divided by b  left associative'#13#10,
  366.   '      a+b - a plus b     a-b - a minus b       left associative'#13#10+
  367.   '    Functions:'#13#10+
  368.   '      Sin(a), Cos(a), Ln(a), Exp(a), Sqrt(a), Atan(a), Abs(a),'#13#10+
  369.   '      Min(a,b), Max(a,b).  Case is ignored.'#13#10+
  370.   'Example:'#13#10+
  371.   '  3df /x1 /y1.2 /s2 min(2,sqrt(abs(1 - x^2-y^2)))'#13#10#10,
  372.   'This is freeware. I ask only for credit and bug reports.'#13#10+
  373.   '    MAJ Gene Ressler  (ressler@cs.cornell.edu)'#13#10+
  374.   '    124 Pine Tree Rd, Ithaca, NY 14580'#13#10);
  375. var
  376.   i: Integer;
  377. begin
  378.   WriteLn(msg+'.');
  379.   for i := 1 to nLines do Write(verbiage[i]);
  380.   Halt;
  381. end;
  382.  
  383. function GetRealParam(name: String; var i: Integer): Real;
  384. var
  385.   x: Real;
  386.   code: Integer;
  387.   param: String;
  388. begin
  389.   param := Copy(ParamStr(i), 3, 255);
  390.   if Length(param) > 0 then
  391.     Val(param, x, code)
  392.   else begin
  393.     inc(i); if i > ParamCount then Usage('Missing '+name);
  394.     Val(ParamStr(i), x, code);
  395.   end;
  396.   if code <> 0 then Usage('Bad '+name);
  397.   GetRealParam := x;
  398. end;
  399.  
  400. procedure GetParams(var expr: String; var cx, cy, size: Real);
  401. var
  402.   i, n, code: Integer;
  403.   param: String[2];
  404. begin
  405.   n := ParamCount;
  406.  
  407.   { Demo case. }
  408.   if n = 0 then begin
  409.     expr := '0.25*cos(28*x^2 + 20*y^2) / (9*x^2 + 6*y^2 + 1)';
  410.     Writeln('/? for help.  Running demo F: ', expr);
  411.     cx := 0.3;
  412.     cy := 0.2;
  413.     size := 1;
  414.     Exit;
  415.   end;
  416.  
  417.   expr := '';
  418.   cx := 0;
  419.   cy := 0;
  420.   size := 1;
  421.  
  422.   i := 1;
  423.   repeat
  424.     param := Copy(ParamStr(i), 1, 2);
  425.     if (param[1] = '/') then
  426.       if Length(param) > 1 then
  427.         case param[2] of
  428.           '?': Usage('Help..');
  429.           's': begin
  430.                  size := GetRealParam('size', i);
  431.                  if size <= 0 then Usage('Non-positive size');
  432.                end;
  433.           'x': cx   := GetRealParam('x', i);
  434.           'y': cy   := GetRealParam('y', i);
  435.           else Usage('Unknown switch '+param);
  436.         end
  437.       else Usage('Missing switch character')
  438.     else begin
  439.       repeat
  440.         if Length(expr) > 0 then expr := expr+' ';
  441.         expr := expr + ParamStr(i);
  442.         inc(i);
  443.       until i > n;
  444.       if expr[1] = '"' then begin
  445.         expr := Copy(expr, 2, 255);
  446.         if expr[Length(expr)] = '"' then
  447.           dec(expr[0])
  448.         else
  449.           Usage('Missing "');
  450.       end;
  451.       Exit;  { Normal exit. }
  452.     end;
  453.     inc(i);
  454.   until i > n;
  455.   Usage('Missing expression');
  456. end;
  457.  
  458. procedure Main;
  459. var
  460.   cx, cy, size: Real;
  461.   expr: String;
  462. begin
  463.   WriteLn('3dF v0.1  MAJ Gene Ressler  April 91');
  464.   GetParams(expr, cx, cy, size);
  465.   ComputeFun(expr, cx, cy, size);
  466.   DrawIt(cx, cy, size);
  467. end;
  468.  
  469. begin
  470.   Main;
  471. end.
  472.