Back to index

radiance  4R0+20100331
gencat.c
Go to the documentation of this file.
00001 #ifndef lint
00002 static const char    RCSid[] = "$Id: gencat.c,v 2.6 2007/09/20 21:28:35 greg Exp $";
00003 #endif
00004 /*****************************************************************************
00005   This program is to make series of right triangles forming hyperbolic cosin
00006   (ie, cosh) curve in between of 2 points.
00007        ^             
00008    f(h)| 
00009        |                  
00010        |               0\ pc:(hc, fc)
00011        |               |    \ 
00012        |               |       \
00013        |             |          \  
00014        |   pa:(ha, fa) 0-------------0 pb:(hb, fb)
00015        |
00016        0------------------------------------------------> h
00017    
00018   Given arguments: 
00019                   material   
00020                   name
00021                   (x0, y0, z0), (x1, y1, z1)
00022                   k        const. value K
00023                   d        distant length desired between 2 points
00024 
00025 ******************************************************************************/
00026 
00027 #include <stdio.h>
00028 #include <stdlib.h>
00029 #include <math.h>
00030 
00031 char  *cmtype, *cname;
00032 double z0, z1;
00033 double k, D;
00034 double d;
00035 double z, h;
00036 
00037 #ifdef notdef
00038 double Newton( b)
00039 double b;
00040 {
00041    if (fabs(cosh(k*D+b)-cosh(b)-(z1-z0)/k) < .001)
00042       return (b);
00043    else {
00044       b = b - (cosh(k*D+b)-cosh(b)-(z1-z0)/k)/(sinh(k*D+b)-sinh(b));
00045       Newton (b);
00046    }
00047 }
00048 #endif
00049 
00050 double Newton(bl)
00051 double bl;
00052 {
00053        double b;
00054        int n = 10000;
00055 
00056        while (n--) {
00057               b = bl- (cosh(k*D+bl)-cosh(bl)-(z1-z0)/k)/(sinh(k*D+bl)-sinh(bl));
00058               if (fabs(b-bl) < .00001)
00059                      return(b);
00060               bl = b;
00061        }
00062 fprintf(stderr, "Interation limit exceeded -- invalid K value\n");
00063        exit(1);
00064 }
00065 
00066 
00067 static void
00068 printhead(ac, av)           /* print command header */
00069 register int  ac;
00070 register char  **av;
00071 {
00072        putchar('#');
00073        while (ac--) {
00074               putchar(' ');
00075               fputs(*av++, stdout);
00076        }
00077        putchar('\n');
00078 }
00079 
00080 int
00081 main (argc, argv)
00082 int argc;
00083 char *argv[];
00084 {
00085    double x0, y0;
00086    double x1, y1;
00087    double b;
00088    double delh;
00089    double f, fprim;
00090    double hb, hc, fb, fc;
00091    int  n;
00092 
00093    if (argc != 11) {
00094       fprintf(stderr, "Usage: %s  material name x0 y0 z0 x1 y1 z1 k d\n",
00095               argv[0]);
00096       exit(1);
00097    }      
00098 
00099    cmtype = argv[1];
00100    cname = argv[2];
00101    x0 = atof(argv[3]); y0 = atof(argv[4]); z0 = atof(argv[5]);
00102    x1 = atof(argv[6]); y1 = atof(argv[7]); z1 = atof(argv[8]);
00103    k = atof(argv[9]); d = atof(argv[10]);
00104    D = sqrt((x1-x0)*(x1-x0) + (y1-y0)*(y1-y0));
00105    b = Newton(0.0); 
00106    z = z0 - k * cosh(b);
00107    printhead(argc, argv);
00108 
00109    n = 0;
00110    for (h=0; h<=D; ) {
00111       f =  k * cosh(k*h+b) + z;
00112       fprim = k* k * sinh(k*h+b);
00113       delh = d / sqrt(1+fprim*fprim); 
00114       fprim = k * k * sinh(k*(h+delh/2)+b); 
00115       hb = sqrt(.01*fprim*fprim/(1+fprim*fprim));
00116       fb = sqrt(.01/(1+(fprim*fprim)));
00117       hc = sqrt(.04/(1+fprim*fprim));
00118       fc = sqrt(.04*fprim*fprim/(1+fprim*fprim));
00119 
00120       printf("\n%s polygon %s.%d\n", cmtype, cname, ++n);
00121       printf("0\n0\n9\n");
00122       printf("%f %f %f\n", h*(x1-x0)/D+x0, h*(y1-y0)/D+y0, f);
00123       if (fprim < 0)
00124       {
00125          printf("%f %f %f\n", (h+hc)*(x1-x0)/D+x0, (h+hc)*(y1-y0)/D+y0, f-fc);
00126          printf("%f %f %f\n", (h+hb)*(x1-x0)/D+x0, (h+hb)*(y1-y0)/D+y0, f+fb);
00127       }
00128       else if (fprim > 0)
00129       {
00130          printf("%f %f %f\n", (h+hc)*(x1-x0)/D+x0, (h+hc)*(y1-y0)/D+y0, f+fc);
00131          printf("%f %f %f\n", (h-hb)*(x1-x0)/D+x0, (h-hb)*(y1-y0)/D+y0, f+fb);
00132       }
00133       else 
00134       { 
00135          printf("%f %f %f\n", (h+.2)*(x1-x0)/D+x0, (h+.2)*(y1-y0)/D+y0, f);
00136          printf("%f %f %f\n", h*(x1-x0)/D+x0, h*(y1-y0)/D+y0, f+.1);
00137       }
00138       h += delh;
00139    }
00140    return 0;
00141 }
00142