home *** CD-ROM | disk | FTP | other *** search
/ rtsi.com / 2014.01.www.rtsi.com.tar / www.rtsi.com / OS9 / OSK / EFFO / forum23.lzh / f23b / SOFTWARE / PDRAW / main.c < prev    next >
C/C++ Source or Header  |  1992-01-15  |  15KB  |  611 lines
  1. /***    main.c   ***/
  2. #define MAIN
  3. #include <stdio.h>
  4. #include <strings.h>
  5. #include <math.h>
  6. #include "header.h"
  7.  
  8. segmptr segmhead;                /* pointer to first segment */
  9. segmptr segmtail;                /* pointer to last segment */
  10.  
  11. pfileptr pfilehead;
  12. pfileptr pfiletail;
  13.  
  14. /* some plot options */
  15. char    *term;                  /* terminal type */
  16. char    xlabel[MAXCHAR];        /* x-label */
  17. char    ylabel[MAXCHAR];        /* y-label */
  18. char    zlabel[MAXCHAR];        /* z-label */
  19. char    toplabel[MAXCHAR];      /* header label */
  20. char    printer[MAXCHAR];       /* printer type */
  21.  
  22. /* info for scaling and transformations */
  23. xyzdata eyepos;                 /* position of the eye/camera */
  24. xyzdata viewcenter;             /* center of view */
  25. xyzdata viewup;                 /* up direction of view */
  26. box     window;                 /* window dimensions relative to center */
  27. box     viewport;               /* viewport in normalized device coordinates*/
  28.  
  29. matrix  I;                      /* identity matrix */
  30. matrix  view_transfo;           /* matrix for view transformation */
  31. double  ndc_width;              /* width in ndc space */
  32. double  ndc_height;             /* height in ndc space */
  33.  
  34. /* plotting options/parameters */
  35. int grid,equalscale,postscript,printplot,noplot;
  36. int xticks,yticks,zticks;
  37. int line,linechange,marker,markerchange;
  38. int linetype[MAXTYPE], markertype[MAXTYPE];
  39. int hiddenline, quick_sort, nosort;
  40. double scale;
  41.  
  42. /* plot boundaries */
  43. double xmin, xmax, ymin, ymax, zmin, zmax;
  44. double pxmin, pymin, pxmax, pymax;
  45.  
  46. main(argc,argv)
  47. int argc;
  48. char *argv[];
  49. {
  50.    FILE     *fopen(), *fp;
  51.    extern   pfileptr pfilehead;
  52.    pfileptr P;
  53.    char     c;
  54.    int      format;
  55.  
  56.    fprintf(stdout,"Pdraw V1.4  9/4/90\n");
  57.    initialize_data();       /* Initialize the variables */
  58.    fileopen(argc,argv);     /* read the data and options */
  59.    initialize_view();       /* initialize the view parameters */
  60.  
  61.    /* process each plotfile */
  62.    for (P=pfilehead; P!=NULL; P=P->next) {
  63.       open_file(P->plotfile,&fp,&format);
  64.  
  65.       fileread(fp,format);     /* read the plotfile */
  66.       read_view(0.0,0.0);      /* read the view transformation */
  67.       draw_view(argc,argv);    /* draw the picture */
  68.       clear_mem();             /* clear the memory */
  69.  
  70.       /* keep on reading until the next blank line */
  71.       while ( (c=getc(fp)) != EOF && c!= '\n');
  72.       while (c!=EOF) {
  73.          if ( (c=getc(fp))!=EOF && c=='\n'){
  74.             /* more curves coming up */
  75.             fileread(fp,format);     /* read the plotfile */
  76.             read_view(0.0,0.0);      /* read the view transformation */
  77.             draw_view(argc,argv);    /* draw the picture */
  78.             clear_mem();             /* clear the memory */
  79.          }
  80.          /* read the rest of the line */
  81.          while ( (c=getc(fp)) != EOF && c!= '\n');
  82.       }
  83.       fclose(fp);
  84.    }
  85.  
  86.    /* read view is used thru every loop iteration to set the view
  87.     * parameters.  draw_view() also calls read_view().  The points 
  88.     * are transformed on a point-by-point basis just before they are
  89.     * plotted using transform_point().  The point-segm structure read
  90.     * in is unaffected by the subsequent matrix operations.
  91.    */
  92.    fprintf(stdout,"Congratulations!!!  You have been TERMINATED...\n");
  93. }
  94.  
  95. /* clear the memory */
  96. clear_mem()
  97. {
  98.    void delete_node();
  99.    void delete_segm();
  100.  
  101.    segmptr S;
  102.    extern  segmptr segmhead,segmtail;
  103.  
  104.    /* delete the nodes in each segment */
  105.    for (S=segmhead; S!=NULL; S=S->next)
  106.       while (S->head!=NULL) delete_node(S->tail,S);
  107.  
  108.    /* now delete all the segments */
  109.    while (segmhead!=NULL) delete_segm(segmtail);
  110.  
  111.    segmhead = NULL;
  112.    segmtail = NULL;
  113. }
  114.  
  115. /* draw the view */
  116. draw_view(argc,argv)
  117. int argc;
  118. char *argv[];
  119. {
  120.    FILE *fopen(), *ips;
  121.    char *sprintf();
  122.    char c, command[MAXCHAR];
  123.    int  PRINT = FALSE;
  124.  
  125.    if (!noplot) {
  126.       if ( (strcmp("xterm",term)==0)  ||
  127.            (strcmp("xterms",term)==0) ||
  128.            (strcmp("vt200",term)==0)  ||
  129.            (strcmp("vt220",term)==0)  ||
  130.            (strcmp("vs100s",term)==0) ) {
  131.          InitWindow(argc,argv);
  132.          RunOps();
  133.       }
  134.    }
  135.  
  136.    /* Prepare for postscript output */
  137.    if (postscript == ON) {
  138.       if ((ips = fopen(PSFILE,"w")) == NULL) {
  139.          fprintf(stderr,"cat: can't open %s\n",PSFILE);
  140.          exit(1);
  141.          }
  142.       initps(ips);
  143.       axesps(ips);
  144.       plotps(ips);
  145.       endgrps(ips);
  146.    }
  147.  
  148.    /* Get postscript printout */
  149.    if (postscript == ON && !printplot) {
  150.       fprintf(stdout,"Type return (or n/q) to ");
  151.       fprintf(stdout,"send the plot to the %s printer :",printer);
  152.       if ( (c=getc(stdin)) == '\n' || c == 'y' || c == 'Y' ) {
  153.          PRINT = ON;
  154.       } else if (c == 'q')
  155.          exit(1);
  156.       if (c != '\n') while ( (c=getc(stdin)) != '\n');
  157.       }
  158.  
  159.    /* now send to the printer */
  160.    if ((printplot || PRINT) && postscript==ON) {
  161.       sprintf(command,"%s %s %s","lpr",printer,PSFILE);
  162.       system(command);
  163.       fprintf(stdout,"Postscript file has been ");
  164.       fprintf(stdout,"sent to the %s printer\n",printer);
  165.    }
  166.  
  167.    if (strcmp("sun",term)==0) {
  168.       fprintf(stdout,"SUN plotting routines not implemented yet.\n");
  169.    }
  170. }
  171.  
  172. /* initialize some data */
  173. initialize_data()
  174. {
  175.    char *getenv();
  176.    char *strcpy();
  177.    char *sprintf();
  178.    int  i,j;
  179.    char *ptr;
  180.  
  181.    /* Identity matrix */
  182.    for (i=0; i<N; i++)
  183.       for (j=0; j<N; j++)
  184.          if (i==j) I[i][j] = 1.0;
  185.          else      I[i][j] = 0.0;
  186.  
  187.    /* Terminal type */
  188.    term = getenv("TERM");
  189.    fprintf(stdout,"   Terminal Type             = %s\n",term);
  190.  
  191.    /* Labeling information */
  192.    strcpy(xlabel,"X-Axis");
  193.    strcpy(ylabel,"Y-Axis");
  194.    strcpy(zlabel,"Z-Axis");
  195.    strcpy(toplabel,"3D Line Plot");
  196.  
  197.    /* plotting parameters/options */
  198.    grid         = OFF;
  199.    equalscale   = ON;
  200.    postscript   = ON;
  201.    printplot    = OFF;
  202.    noplot       = OFF;
  203.    line         = ON;
  204.    linechange   = OFF;
  205.    marker       = OFF;
  206.    markerchange = OFF;
  207.    hiddenline   = OFF;
  208.    quick_sort   = ON;
  209.    nosort       = OFF;
  210.    for (i=0; i<MAXTYPE; i++) linetype[i] = i;
  211.    for (i=0; i<MAXTYPE; i++) markertype[i] = i;
  212.    scale        = 1.0;
  213.    xticks       = 2;
  214.    yticks       = 2;
  215.    zticks       = 2;
  216.    
  217.    /* data */
  218.    segmhead     = NULL;
  219.    segmtail     = NULL;
  220.  
  221.    /* printer type */
  222.    if ((ptr = getenv("PRINTER")) == NULL)
  223.       strcpy(printer,"-Plp550M");
  224.    else
  225.       sprintf(printer,"-P%s",ptr);
  226.    fprintf(stdout,"   Default Printer           = %s\n",printer);
  227.  
  228.    /* ndc parameters */
  229.    ndc_width = 1.0;
  230.    ndc_height = 1.0;
  231. }
  232.  
  233. initialize_view()
  234. {
  235.    extern double xmax, xmin, ymax, ymin, zmax, zmin;
  236.    double max;
  237.  
  238.    /* viewing parameters */
  239.    if ((eyepos.x*eyepos.x + eyepos.y*eyepos.y + eyepos.z*eyepos.z) < SMALL) {
  240.       eyepos.x  = 1.0;
  241.       eyepos.y  = 1.5;
  242.       eyepos.z  = 0.5;
  243.    }
  244.  
  245.    viewcenter.x  = 0.0;
  246.    viewcenter.y  = 0.0;
  247.    viewcenter.z  = 0.0;
  248.  
  249.    viewup.x = 0.0;
  250.    viewup.y = 0.0;
  251.    viewup.z = 1.0;
  252.  
  253.    max = -1.0e10;
  254.    if (xmax-xmin > max) max = xmax-xmin;
  255.    if (ymax-ymin > max) max = zmax-ymin;
  256.    if (zmax-zmin > max) max = zmax-zmin;
  257.    max = 1.3*max;
  258.    window.xl = -max;
  259.    window.xr =  max;
  260.    window.yb = -max;
  261.    window.yt =  max;
  262.  
  263.    viewport.xl = 1.0*PSMIN;
  264.    viewport.xr = 1.0*PSMAX*scale;
  265.    viewport.yb = 1.0*PSMIN;
  266.    viewport.yt = 1.0*PSMAX*scale;
  267. }
  268.  
  269. read_view(theta,phi)
  270. double theta, phi;
  271. {
  272.    xyzdata transform_point(), normalize_vector();
  273.    matrix vwscale, R, Rp, Rx, Ry, Rz, SC;
  274.    double dxsq, dysq, dzsq, l, r, v, angle, theta0, phi0;
  275.    xyzdata newpoint, neweye;
  276.  
  277.    /* eye/center */
  278.    dxsq = (eyepos.x-viewcenter.x)*(eyepos.x-viewcenter.x);
  279.    dysq = (eyepos.y-viewcenter.y)*(eyepos.y-viewcenter.y);
  280.    dzsq = (eyepos.z-viewcenter.z)*(eyepos.z-viewcenter.z);
  281.    l = sqrt(dxsq + dysq + dzsq);
  282.    r = sqrt(dxsq + dysq);
  283.    if (l < SMALL) {
  284.       fprintf(stderr,"Error: Zero Length Vector - No View Specified\n");
  285.       exit(1);
  286.    }
  287.    phi0   = acos((eyepos.z-viewcenter.z)/l);
  288.    if (r > SMALL) theta0 = acos((eyepos.x-viewcenter.x)/r);
  289.    else           theta0 = 0.0;
  290.    neweye.x = l*sin(phi0+phi)*cos(theta0+theta);
  291.    neweye.y = l*sin(phi0+phi)*sin(theta0+theta);
  292.    neweye.z = l*cos(phi0+phi);
  293.  
  294.    /* eye/center */
  295.    dxsq = (neweye.x-viewcenter.x)*(neweye.x-viewcenter.x);
  296.    dysq = (neweye.y-viewcenter.y)*(neweye.y-viewcenter.y);
  297.    dzsq = (neweye.z-viewcenter.z)*(neweye.z-viewcenter.z);
  298.    l = sqrt(dxsq + dysq + dzsq);
  299.    v = sqrt(       dysq + dzsq);
  300.    if (l==0) {
  301.       fprintf(stderr,"Error: Zero Length Vector - No View Specified\n");
  302.       exit(1);
  303.    }
  304.  
  305.    /* Rotate around x by Theta */
  306.    if (v > SMALL) angle = acos((neweye.z-viewcenter.z)/v); 
  307.    else           angle = 0.0;
  308.    if (neweye.y < 0.0) angle = -1.0*angle;
  309.    rotate('x',angle,Rx);
  310.  
  311.    /* Rotate around y by -Theta */
  312.    angle = -acos(v/l); 
  313.    if (neweye.x < 0.0) angle = -1.0*angle;
  314.    rotate('y',angle,Ry);
  315.  
  316.    /* Rotation Matrix */
  317.    mult_matrix(Rp,Rx,Ry);                     /* Rp = Rx * Ry */
  318.  
  319.    /* Find the necessary matrix to rotate to the view up position */
  320.    if ( dxsq + dysq > SMALL) {
  321.       newpoint   = transform_point(viewup,Rp);
  322.       newpoint.z = 0.0;
  323.       newpoint   = normalize_vector(newpoint);
  324.       angle      = -acos(newpoint.y);  /* dot_product[(newpoint),(0,1,0)] */
  325.       if (newpoint.x>0.0)
  326.          angle = -1.0*angle;           /* from cross_product with (0,1,0) */
  327.    } else {
  328.       angle = 0.0;
  329.    }
  330.    rotate('z',angle,Rz);
  331.  
  332.    /* Rotation Matrix */
  333.    mult_matrix(R,Rp,Rz);                     /* R = R * Ry */
  334.  
  335.    /* xyz scaling */
  336.    if (equalscale==OFF) {
  337.       assign_matrix(SC,I);
  338.       if ((xmax-xmin) > 1.0e-10) SC[0][0] = 1.0/(xmax-xmin);
  339.       if ((ymax-ymin) > 1.0e-10) SC[1][1] = 1.0/(ymax-ymin);
  340.       if ((zmax-zmin) > 1.0e-10) SC[2][2] = 1.0/(zmax-zmin);
  341.       /* multiply the rotation matrix by a scaling matrix */
  342.       mult_matrix(R,SC,R);
  343.    }
  344.  
  345.    /* window - viewport scaling */
  346.    assign_matrix(vwscale,I);
  347.  
  348.    find_window(R);
  349.    vwscale[0][0] = (viewport.xr - viewport.xl)/(window.xr-window.xl);
  350.    vwscale[1][1] = (viewport.yt - viewport.yb)/(window.yt-window.yb);
  351.    vwscale[3][0] = viewport.xl - vwscale[0][0]*window.xl;
  352.    vwscale[3][1] = viewport.yb - vwscale[1][1]*window.yb;
  353.  
  354.    /* Apply scaling and transformation */
  355.    mult_matrix(view_transfo,R,vwscale);
  356. }
  357.  
  358. find_window(R)
  359. matrix R;
  360. {
  361.    xyzdata transform_point();
  362.    xyzdata point;
  363.    double  delx, dely;
  364.    int     i,j,k;
  365.  
  366.    window.xl =  1.0e10;
  367.    window.xr = -1.0e10;
  368.    window.yb =  1.0e10;
  369.    window.yt = -1.0e10;
  370.  
  371.    for (i=0; i<2; i++)
  372.    for (j=0; j<2; j++)
  373.    for (k=0; k<2; k++) {
  374.       point.x = xmin + i*(xmax-xmin);
  375.       point.y = ymin + j*(ymax-ymin);
  376.       point.z = zmin + k*(zmax-zmin);
  377.       point = transform_point(point,R);
  378.       if (point.x < window.xl) window.xl = point.x;
  379.       if (point.x > window.xr) window.xr = point.x;
  380.       if (point.y < window.yb) window.yb = point.y;
  381.       if (point.y > window.yt) window.yt = point.y;
  382.    }
  383.  
  384.    delx = window.xr - window.xl;
  385.    dely = window.yt - window.yb;
  386.    if (delx > dely) window.yt = window.yb + delx;
  387.    if (dely > delx){
  388.       window.xl = window.xl - 0.5*(dely-delx);
  389.       window.xr = window.xr + 0.5*(dely-delx);
  390.    }
  391. }
  392.  
  393. apply_view()
  394. {
  395.    extern segmptr segmhead;
  396.  
  397.    double a[N];
  398.    int    j;
  399.    segmptr S;
  400.    nodeptr Nd;
  401.  
  402.    for (S=segmhead; S!=NULL; S=S->next)
  403.       for (Nd=S->head; Nd!=NULL; Nd=Nd->next) {
  404.          for (j=0; j<N; j++)
  405.             a[j] = Nd->x * view_transfo[0][j] + 
  406.                    Nd->y * view_transfo[1][j] + 
  407.                    Nd->z * view_transfo[2][j] +
  408.                    view_transfo[3][j];
  409.          Nd->x = a[0]/a[3];
  410.          Nd->y = a[1]/a[3];
  411.          Nd->z = a[2]/a[3];
  412.       }
  413. }
  414.    
  415. check_view_angles(theta0,phi0)
  416. double *theta0, *phi0;
  417. {
  418.    extern xyzdata eyepos, viewcenter;
  419.  
  420.    double dx, dy, dz, l, r;
  421.  
  422.    /* view vector */
  423.    dx = eyepos.x - viewcenter.x;
  424.    dy = eyepos.y - viewcenter.y;
  425.    dz = eyepos.z - viewcenter.z;
  426.    l  = sqrt(dx*dx + dy*dy + dz*dz);
  427.    r  = sqrt(dx*dx + dy*dy);
  428.    if (l < SMALL) {
  429.       fprintf(stderr,"Error: Zero Length Vector - No View Specified\n");
  430.       exit(1);
  431.    }
  432.  
  433.    /* angle from z-axis */
  434.    *phi0 = acos(dz/l);
  435.  
  436.    /* angle from x-axis */
  437.    *theta0 = 0.0;
  438.    if (r > SMALL) *theta0 = acos(dx/r);
  439. }
  440.  
  441. xyzdata transform_point(point,A)
  442. xyzdata point; 
  443. matrix A;
  444. {
  445.    double a[N];
  446.    int    j;
  447.    xyzdata newpoint;
  448.  
  449.    for (j=0; j<N; j++)
  450.       a[j] = point.x * A[0][j] + 
  451.              point.y * A[1][j] + 
  452.              point.z * A[2][j] + A[3][j];
  453.    newpoint.x = a[0]/a[3];
  454.    newpoint.y = a[1]/a[3];
  455.    newpoint.z = a[2]/a[3];
  456.    return(newpoint);
  457. }
  458.  
  459. xyzdata normalize_vector(vector)
  460. xyzdata vector;
  461. {
  462.       double dl;
  463.       xyzdata newvector;
  464.  
  465.       dl = sqrt(vector.x*vector.x + vector.y*vector.y + vector.z*vector.z);
  466.       newvector.x = vector.x/dl;
  467.       newvector.y = vector.y/dl;
  468.       newvector.z = vector.z/dl;
  469.       return(newvector);
  470. }
  471.  
  472. xyzdata midpoint(x1,y1,z1,x2,y2,z2,ratio)
  473. double  x1,y1,z1,x2,y2,z2,ratio;
  474. {
  475.    xyzdata midpt, newpt;
  476.    xyzdata transform_point();
  477.  
  478.    midpt.x = x1+ratio*(x2-x1);
  479.    midpt.y = y1+ratio*(y2-y1);
  480.    midpt.z = z1+ratio*(z2-z1);
  481.  
  482.    newpt = transform_point(midpt,view_transfo);
  483.    return(newpt);
  484. }
  485.  
  486. print_point(point)
  487. xyzdata point;
  488. {
  489.       fprintf(stdout,"x = %10.5f ", point.x);
  490.       fprintf(stdout,"y = %10.5f ", point.y);
  491.       fprintf(stdout,"z = %10.5f\n",point.z);
  492. }
  493.  
  494. print_data()
  495. {
  496.    extern segmptr segmhead;
  497.  
  498.    segmptr S;
  499.    nodeptr Nd;
  500.  
  501.    for (S=segmhead; S!=NULL; S=S->next)
  502.       for (Nd=S->head; Nd!=NULL; Nd=Nd->next) {
  503.          fprintf(stdout,"x = %10.5f ",Nd->x);
  504.          fprintf(stdout,"y = %10.5f ",Nd->y);
  505.          fprintf(stdout,"z = %10.5f\n",Nd->z);
  506.       }
  507. }
  508.  
  509. /* matrix OPERATIONS */
  510. print_matrix(A)
  511. matrix A;
  512. {
  513.    int i,j;
  514.  
  515.    for (i=0; i<N; i++) {
  516.       fprintf(stdout,"\n");
  517.       for (j=0; j<N; j++) 
  518.          fprintf(stdout,"%10.5f",A[i][j]);
  519.    }
  520.    fprintf(stdout,"\n");
  521. }
  522.  
  523. assign_matrix(A,B)
  524. matrix A,B;
  525. {
  526.    int i,j;
  527.  
  528.    for (i=0; i<N; i++)
  529.    for (j=0; j<N; j++)
  530.       A[i][j] = B[i][j];
  531. }
  532.  
  533. mult_matrix(A,B,C)
  534. matrix A,B,C;
  535. {
  536.    int i,j;
  537.  
  538.    for (i=0; i<N; i++)
  539.    for (j=0; j<N; j++)
  540.       A[i][j] = B[i][0]*C[0][j] + B[i][1]*C[1][j] +
  541.                 B[i][2]*C[2][j] + B[i][3]*C[3][j];
  542. }
  543.  
  544. scalematrix(sx,sy,sz,A)
  545. double sx, sy, sz;
  546. matrix A;
  547. {
  548.    assign_matrix(A,I);
  549.    A[0][0] = sx;
  550.    A[1][1] = sy;
  551.    A[2][2] = sz;
  552. }
  553.  
  554. mirror(axis,A)
  555. char axis;
  556. matrix A;
  557. {
  558.    double sx, sy, sz;
  559.  
  560.    sx = 1.0;
  561.    sy = 1.0;
  562.    sz = 1.0;
  563.    if (axis == 'x') sx = -1.0;
  564.    else if (axis == 'y') sy = -1.0;
  565.    else if (axis == 'z') sz = -1.0;
  566.    else {
  567.     fprintf(stderr,"Error: The character %c does not belong in mirror\n",axis);
  568.     exit(1);
  569.    }
  570.  
  571.    scalematrix(sx,sy,sz,A);
  572. }
  573.  
  574. rotate(axis,angle,A)
  575. char axis;
  576. double angle;
  577. matrix A;
  578. {   
  579.    assign_matrix(A,I);
  580.    if (axis == 'z') {
  581.       A[0][0] =  cos(angle);
  582.       A[0][1] =  sin(angle);
  583.       A[1][0] = -sin(angle);
  584.       A[1][1] =  cos(angle);
  585.    } else if (axis == 'y') {
  586.       A[0][0] =  cos(angle);
  587.       A[0][2] = -sin(angle);
  588.       A[2][0] =  sin(angle);
  589.       A[2][2] =  cos(angle);
  590.    } else if (axis == 'x') {
  591.       A[1][1] =  cos(angle);
  592.       A[1][2] =  sin(angle);
  593.       A[2][1] = -sin(angle);
  594.       A[2][2] =  cos(angle);
  595.    } else {
  596.     fprintf(stderr,"Error: The character %c does not belong in rotate\n",axis);
  597.     exit(1);
  598.    }
  599. }
  600.  
  601. translate(tx,ty,tz,A)
  602. double tx,ty,tz;
  603. matrix A;
  604. {
  605.    assign_matrix(A,I);
  606.    A[3][0] = tx;
  607.    A[3][1] = ty;
  608.    A[3][2] = tz;
  609. }
  610.  
  611.