home *** CD-ROM | disk | FTP | other *** search
/ Otherware / Otherware_1_SB_Development.iso / mac / developm / source / graph3d.sit / g3D_CalcFunc.c < prev    next >
Encoding:
C/C++ Source or Header  |  1989-08-28  |  10.9 KB  |  605 lines
  1. /*
  2.     Copyright '89    Christopher Moll
  3.     all rights reserved
  4. */
  5.  
  6.  
  7. #include    "graph3D.h"
  8. #include    <math.h>
  9.  
  10. #ifndef    _LSC3_
  11. #    include    <QuickDraw.h>
  12. #    include    <MemoryMgr.h>
  13. #endif
  14.  
  15.  
  16. extern    Boolean        useMainFunc;
  17.  
  18. extern    int        mainFunct[64];
  19. extern    int        numMnOps;
  20. extern    Real    mnConsts[30];
  21.  
  22. extern    int        derivFunct[64];
  23. extern    int        numDerivOps;
  24. extern    Real    derivConsts[30];
  25.  
  26.  
  27. Real
  28. TheFunction(x, y, z)
  29. Real    x, y, z;
  30. {
  31.     Real    stack[50];
  32. register    int        i, currOp;
  33.     Real        StackOp(), FuncEval(), ArithEval();
  34. register    Real    *stackPtr;
  35.     int        currConst = 0;
  36.     int        numOps;
  37.  
  38.     stackPtr = stack + 49;
  39.     if (useMainFunc)
  40.         if (numMnOps EQ 0)
  41.             return(0.0);
  42.         else
  43.             numOps = numMnOps;
  44.     else
  45.         if (numDerivOps EQ 0)
  46.             return(0.0);
  47.         else
  48.             numOps = numDerivOps;
  49.  
  50.     *stackPtr = 0.0;
  51.     for (i = 0; i < numOps; i++)
  52.     {
  53.         if (useMainFunc)
  54.             currOp = mainFunct[i];
  55.         else
  56.             currOp = derivFunct[i];
  57.         if (IS_PUSH(currOp))
  58.         {
  59.             --stackPtr;
  60.             switch (currOp) {
  61.             case PUSH_X_OP:
  62.                 *stackPtr = x;
  63.                 break;
  64.             case PUSH_Y_OP:
  65.                 *stackPtr = y;
  66.                 break;
  67.             case PUSH_Z_OP:
  68.                 *stackPtr = z;
  69.                 break;
  70.             case PUSH_NUM_OP:
  71.                 if (useMainFunc)
  72.                     *stackPtr = mnConsts[currConst++];
  73.                 else
  74.                     *stackPtr = derivConsts[currConst++];
  75.                 break;
  76.             }
  77.         }
  78.         else if (IS_FUNC(currOp))
  79.             *stackPtr = FuncEval(*stackPtr, currOp);
  80.         else
  81.         {
  82.             stackPtr[1] = ArithEval(stackPtr[1], stackPtr[0], currOp);
  83.             ++stackPtr;
  84.         }
  85.     }
  86.  
  87.     return(*stackPtr);
  88. }
  89.  
  90. static
  91. Real
  92. ArithEval(arg1, arg2, op)
  93. Real    arg1, arg2;
  94. int        op;
  95. {
  96. extern    int    graphType;
  97.  
  98.     switch (op) {
  99.     case MULT_OP:
  100.         return(arg1 * arg2);
  101.     case DIV_OP:
  102.         if (graphType EQ DERIV_GTYPE)
  103.             if (Abs(arg2) > 1e-50)
  104.                 return(arg1 / arg2);
  105.             else
  106.                 return(0.0);
  107.         else
  108.             return(arg1 / arg2);
  109.     case ADD_OP:
  110.         return(arg1 + arg2);
  111.     case SUBT_OP:
  112.         return(arg1 - arg2);
  113.     case MOD_OP:
  114.         return(round(arg1) MOD round(arg2));
  115.     case POWER_OP:
  116.         return(pow(arg1, arg2));
  117.     }
  118. }
  119.  
  120. kw()
  121. {
  122.     if (NOT(OptionDown()))
  123.     {
  124.         while(NOT(ShiftDown()));
  125.         while(ShiftDown());
  126.     }
  127. }
  128.  
  129. es()
  130. {
  131.     ExitToShell();
  132. }
  133.  
  134. static
  135. Real
  136. FuncEval(arg, op)
  137. Real    arg;
  138. int        op;
  139. {
  140.     switch (op) {
  141.     case NEG_OP:
  142.         return(-arg);
  143.     case EXP_OP:
  144.         return(exp(arg));
  145.     case LOG_OP:
  146.         return(log(arg));
  147.     case LOG10_OP:
  148.         return(log10(arg));
  149.     case SQ_OP:
  150.         return(SQ(arg));
  151.     case SQRT_OP:
  152.         if (arg > 0)
  153.             return(sqrt(arg));
  154.         else
  155.             return(0.0);
  156.     case SIN_OP:
  157.         return(sin(arg));
  158.     case ASIN_OP:
  159.         return(asin(arg));
  160.     case COS_OP:
  161.         return(cos(arg));
  162.     case ACOS_OP:
  163.         return(acos(arg));
  164.     case TAN_OP:
  165.         return(tan(arg));
  166.     case ATAN_OP:
  167.         return(atan(arg));
  168.     case COTAN_OP:
  169.         return(cotan(arg));
  170.     case SINH_OP:
  171.         return(sinh(arg));
  172.     case COSH_OP:
  173.         return(cosh(arg));
  174.     case TANH_OP:
  175.         return(tanh(arg));
  176.     }
  177. }
  178.  
  179. /*-------------------------------------------------------------------------*/
  180.  
  181.  
  182. extern    Real        *funcResults;
  183. extern    Vector        *scalVectResults;
  184. extern    Point        *graphPoints;    /* same elements as funcResults */
  185.  
  186.  
  187. extern    Real    startX, startY;
  188. extern    Real    endX, endY;
  189. extern    Real    deltaX, deltaY;
  190. extern    int        numX, numY;
  191.  
  192. extern    Real    scale;
  193.  
  194. extern    Boolean        vectrsCurrent;
  195. extern    Boolean        pntsCurrent;
  196.  
  197. extern    Vector        maxVect, minVect;
  198. extern    Point    XaxisPt, YaxisPt, ZaxisPt;
  199. extern    Point    Origin;
  200. extern    Real    rotMatrx[3][3];
  201.  
  202. extern    Rect    graphRect;
  203. extern    Boolean        hideSurface, useHeight;
  204.  
  205. extern    Real    zScale;
  206. extern    Real    ceilgZ, floorZ;
  207.  
  208.  
  209. SRotatePoints(cAng, sAng, plane)
  210. register    Real    cAng, sAng;
  211. int        plane;
  212. {
  213.     register    long    i, numElmnts;
  214.     register    double    temp1, temp2;
  215.  
  216.     numElmnts = numX * numY;
  217.     if (plane EQ 1)
  218.     {
  219.         for (i = 0; i < numElmnts; i++)
  220.         {
  221.             temp1 = scalVectResults[i].y;
  222.             temp2 = scalVectResults[i].z;
  223.             scalVectResults[i].y = temp1 * cAng + temp2 * sAng;
  224.             scalVectResults[i].z = temp2 * cAng - temp1 * sAng;
  225.         }
  226.     }
  227.     else if (plane EQ 2)
  228.     {
  229.         for (i = 0; i < numElmnts; i++)
  230.         {
  231.             temp1 = scalVectResults[i].x;
  232.             temp2 = scalVectResults[i].z;
  233.             scalVectResults[i].x = temp1 * cAng - temp2 * sAng;
  234.             scalVectResults[i].z = temp2 * cAng + temp1 * sAng;
  235.         }
  236.     }
  237.     else if (plane EQ 3)
  238.     {
  239.         for (i = 0; i < numElmnts; i++)
  240.         {
  241.             temp1 = scalVectResults[i].x;
  242.             temp2 = scalVectResults[i].y;
  243.             scalVectResults[i].x = temp1 * cAng - temp2 * sAng;
  244.             scalVectResults[i].y = temp2 * cAng + temp1 * sAng;
  245.         }
  246.     }
  247.     else
  248.     {
  249.         for (i = 0; i < numElmnts; i++)
  250.         {
  251.             temp1 = scalVectResults[i].z;
  252.             temp2 = scalVectResults[i].y;
  253.             scalVectResults[i].z = temp1 * cAng - temp2 * sAng;
  254.             scalVectResults[i].y = temp2 * cAng + temp1 * sAng;
  255.         }
  256.     }
  257. }
  258.  
  259. SetPoints()
  260. {
  261. register    Real    _scale;
  262. register    int    xCnt, yCnt;
  263. register    long    offSet;
  264.     
  265.     SetNewScale();
  266.  
  267.     _scale = scale;
  268.     if (hideSurface)
  269.         CheckHideView();
  270.  
  271.     for (xCnt = 0; xCnt < numX; xCnt++)
  272.         for (yCnt = 0; yCnt < numY; yCnt++)
  273.         {
  274.             offSet = (long)xCnt * (long)numY + yCnt;
  275.             graphPoints[offSet].h = (int)(scalVectResults[offSet].y * _scale)
  276.                             + Origin.h;
  277.             graphPoints[offSet].v = (int)(-scalVectResults[offSet].z * _scale)
  278.                             + Origin.v;
  279.         }
  280.     GenAxes();
  281.     pntsCurrent = TRUE;
  282.     if (hideSurface AND useHeight)
  283.         SetHeightShades();
  284.     else if (hideSurface)
  285.         SetNormalShades();
  286. }
  287.  
  288. /*************************/
  289.  
  290. Uchar        *gridFace = NIL;
  291. extern        Real        cosSrc, sinSrc;
  292. extern        Boolean        scrnColor, colorNorms;
  293.  
  294.  
  295. SetNormalShades()
  296. {
  297. register    int        xCnt, yCnt, v;
  298. register    long    offSet;
  299.     Uchar    VectFace();
  300.  
  301.     colorNorms = scrnColor;
  302.     if (gridFace NEQ NIL)
  303.         free(gridFace);
  304.     gridFace = (Uchar *)NewPtr(sizeof(char) * (long)numX * numY);
  305.     if (MemErr)
  306.         return(FALSE);
  307.  
  308.     for (xCnt = 0; xCnt < numX-1; xCnt++)
  309.     {
  310.         for (yCnt = 0; yCnt < numY-1; yCnt++)
  311.         {
  312.             offSet = (long)xCnt * (long)numY + (long)yCnt;
  313.             v = VectFace(&scalVectResults[(long)(xCnt + 1) * numY + (long)yCnt],
  314.                     &scalVectResults[offSet],
  315.                     &scalVectResults[(long)xCnt * numY + (long)(yCnt + 1)]);
  316.             gridFace[offSet] = v;
  317.         }
  318.     }
  319. }
  320.  
  321. Uchar
  322. VectFace(a, b, c)
  323. Vector    *a, *b, *c;
  324. {
  325.     Vector    side1, side2, reslt;
  326.     Real    length, zCoor;
  327.     Real    VectMag();
  328. register    Real    temp1, temp2;
  329.  
  330.     SubtVector(a, b, &side1);
  331.     SubtVector(c, b, &side2);
  332.     CrossProduct(&side1, &side2, &reslt);
  333.  
  334.             temp1 = reslt.y;
  335.             temp2 = reslt.z;
  336.             reslt.y = temp1 * cosSrc + temp2 * sinSrc;
  337.             reslt.z = temp2 * cosSrc - temp1 * sinSrc;
  338.  
  339. /*if (reslt.x < 0)
  340.     return(0);
  341. */
  342.     length = VectMag(&reslt);
  343.     zCoor = reslt.z/length;
  344.     if (zCoor < 0)
  345.         zCoor = 0;
  346.     if (scrnColor)
  347.         zCoor = zCoor * 255.0;
  348.     else
  349.         zCoor = zCoor * 63.0;
  350.  
  351.     return((Uchar)round(zCoor));
  352. }
  353.  
  354.  
  355. /* cross-product of two vectors */
  356. CrossProduct(a, b, reslt)
  357. register    Vector    *a, *b, *reslt;
  358. {
  359.     reslt->x = a->y * b->z - a->z * b->y;
  360.     reslt->y = a->z * b->x - a->x * b->z;
  361.     reslt->z = a->x * b->y - a->y * b->x;
  362. }
  363.  
  364. /*** Difference and sum of two vectors;
  365.             result vector can be same as a or b ***/
  366. SubtVector(a, b, result)
  367. register    Vector    *a, *b, *result;
  368. {
  369.     result->x = a->x - b->x;
  370.     result->y = a->y - b->y;
  371.     result->z = a->z - b->z;
  372. }
  373.  
  374. Real
  375. VectMag(vectP)
  376. register    Vector    *vectP;
  377. {
  378.     return( sqrt( SQ(vectP->x) + SQ(vectP->y) + SQ(vectP->z) ) );
  379. }
  380.  
  381.  
  382. SetHeightShades()
  383. {
  384. register    int        xCnt, yCnt, v;
  385. register    long    offSet;
  386.     Real    min, max, sc;
  387.  
  388.     colorNorms = scrnColor;
  389.     if (gridFace NEQ NIL)
  390.         free(gridFace);
  391.     gridFace = (Uchar *)NewPtr(sizeof(char) * (long)numX * numY);
  392.     if (MemErr)
  393.         return(FALSE);
  394.  
  395.     FindFuncMaxMin(funcResults, &min, &max);
  396.     if (scrnColor)
  397.         sc = 255.0 / (max - min);
  398.     else
  399.         sc = 63.0 / (max - min);
  400.  
  401.     for (xCnt = 0; xCnt < numX-1; xCnt++)
  402.     {
  403.         for (yCnt = 0; yCnt < numY-1; yCnt++)
  404.         {
  405.             offSet = (long)xCnt * (long)numY + (long)yCnt;
  406.             gridFace[offSet] = (funcResults[offSet] - min) * sc;
  407.         }
  408.     }
  409. }
  410.  
  411. static
  412. FindFuncMaxMin(funcList, minPtr, maxPtr)
  413. register    Real    *funcList;
  414. Real        *minPtr, *maxPtr;
  415. {
  416. register    long    i, numElems;
  417. register    Real    max, min, funcVal;
  418.  
  419.     max = min = funcList[0];
  420.  
  421.     numElems = (long)numX * (long)numY;
  422.     for (i = 1; i < numElems; i++)
  423.     {
  424.         funcVal = funcList[i];
  425.         if (funcVal > max)
  426.             max = funcVal;
  427.         else if (funcVal < min)
  428.             min = funcVal;
  429.     }
  430.  
  431.     *minPtr = min;
  432.     *maxPtr = max;
  433. }
  434.  
  435.  
  436. /*************************/
  437.  
  438. SetVectrs()
  439. {
  440. register    int    xCnt, yCnt;
  441. register    long    offSet;
  442.     Real        ImposLimits();
  443. register    Real        xValue, yValue;
  444. register    Vector        *scalePtr;
  445.     Real    SdeltaX;
  446. register    Real    SdeltaY;
  447.     Vector    pntVect;
  448.  
  449.     SdeltaX = deltaX;
  450.     SdeltaY = deltaY;
  451.  
  452.     scalePtr = scalVectResults;
  453.  
  454.     xValue = startX;
  455.  
  456.     for (xCnt = 0; xCnt < numX; xCnt++)
  457.     {
  458.         yValue = startY;
  459.         pntVect.x = xValue;
  460.         for (yCnt = 0; yCnt < numY; yCnt++)
  461.         {
  462.             offSet = (long)xCnt * (long)numY + yCnt;
  463.             pntVect.y = yValue;
  464.             pntVect.z = (ImposLimits(funcResults[offSet]));
  465.             MatrxMul(rotMatrx, &pntVect, scalePtr++);
  466.             yValue += SdeltaY;
  467.         }
  468.         xValue += SdeltaX;
  469.     }
  470.     vectrsCurrent = TRUE;
  471. }
  472.  
  473.  
  474. static
  475. SetNewScale()
  476. {
  477.     DoubPoint    minPnt, maxPnt;
  478.     Real        scaleX, scaleY;
  479.     int            graphWidth, graphHeight;
  480.     DoubPoint    doubOrigin;
  481.  
  482.     FindMaxMin(scalVectResults, &minPnt, &maxPnt);
  483.  
  484.     graphWidth = graphRect.right - graphRect.left;
  485.     graphHeight = graphRect.bottom - graphRect.top;
  486.  
  487.     scaleX = (Real)(graphWidth) / (maxPnt.h - minPnt.h);
  488.     scaleY = (Real)(graphHeight) / (maxPnt.v - minPnt.v);
  489.     scale = MIN(scaleX, scaleY);
  490.  
  491.     maxPnt.h = maxPnt.h * scale;
  492.     maxPnt.v = maxPnt.v * scale;
  493.     minPnt.h = minPnt.h * scale;
  494.     minPnt.v = minPnt.v * scale;
  495.  
  496.     doubOrigin.h = ((graphWidth >> 1) + graphRect.left)
  497.                 - ( ((maxPnt.h - minPnt.h) / 2) + minPnt.h );
  498.     doubOrigin.v = ((graphHeight >> 1) + graphRect.top)
  499.             - ( ((maxPnt.v - minPnt.v) / 2) + minPnt.v );
  500.     Origin.h = round(doubOrigin.h);
  501.     Origin.v = round(doubOrigin.v);
  502. }
  503.  
  504. static
  505. FindMaxMin(graphVects, minPnt, maxPnt)
  506. register    Vector    *graphVects;
  507. DoubPoint    *minPnt, *maxPnt;
  508. {
  509. register    long    i, numElems;
  510. register    Real    max_h, max_v;
  511. register    Real    min_h, min_v;
  512.     Real    vertVal;
  513.  
  514.     max_h = min_h = graphVects[0].y;
  515.     max_v = min_v = -graphVects[0].z;
  516.  
  517.     numElems = (long)numX * (long)numY;
  518.     for (i = 1; i < numElems; i++)
  519.     {
  520.         if (graphVects[i].y > max_h)
  521.             max_h = graphVects[i].y;
  522.         else if (graphVects[i].y < min_h)
  523.             min_h = graphVects[i].y;
  524.         vertVal = -graphVects[i].z;
  525.         if (vertVal > max_v)
  526.             max_v = vertVal;
  527.         else if (vertVal < min_v)
  528.             min_v = vertVal;
  529.     }
  530.  
  531.     minPnt->h = min_h;
  532.     minPnt->v = min_v;
  533.     maxPnt->h = max_h;
  534.     maxPnt->v = max_v;
  535. }
  536.  
  537. Real
  538. ImposLimits(funcVal)
  539. Real    funcVal;
  540. {
  541.     Real    temp;
  542.     
  543.     temp = funcVal;
  544.     
  545.     if (temp > ceilgZ)
  546.         temp = ceilgZ;
  547.     else if (temp < floorZ)
  548.         temp = floorZ;
  549.  
  550.     return(temp * zScale);
  551. }
  552.  
  553. static
  554. ExpandRect(theRect, thePoint)
  555. register    DoubRect    *theRect;
  556. register    DoubPoint    *thePoint;
  557. {
  558.     if (thePoint->h < theRect->left)
  559.         theRect->left = thePoint->h;
  560.     else if (thePoint->h > theRect->right)
  561.         theRect->right = thePoint->h;
  562.  
  563.     if (thePoint->v < theRect->top)
  564.         theRect->top = thePoint->v;
  565.     else if (thePoint->v > theRect->bottom)
  566.         theRect->bottom = thePoint->v;
  567. }
  568.  
  569.  
  570. #ifdef    _MC68881_    /*** This function will be a macro calling ldexp() otherwise ***/
  571.  
  572. Real
  573. mldexp(a, n)
  574. Real    a;
  575. register    int        n;
  576. {
  577.     if (n EQ 0)
  578.         return(a);
  579.     else if (n > 0)
  580.         return(a * (Real)((long)2 << (n - 1)));
  581.     else
  582.         return(a / (Real)((long)2 << -(n + 1)));
  583. }
  584.  
  585. #endif
  586.  
  587. static
  588. GenAxes()
  589. {
  590.     Vector        axisVect;
  591.  
  592.     axisVect.x = 1000;
  593.     axisVect.y = 0;
  594.     axisVect.z = 0;
  595.     Conv3VtoPt(&axisVect, &XaxisPt);
  596.  
  597.     axisVect.x = 0;
  598.     axisVect.y = 1000;
  599.     Conv3VtoPt(&axisVect, &YaxisPt);
  600.  
  601.     axisVect.y = 0;
  602.     axisVect.z = 1000;
  603.     Conv3VtoPt(&axisVect, &ZaxisPt);
  604. }
  605.