ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/SHAPES/visualizer.cpp
Revision: 1290
Committed: Wed Jun 23 20:53:17 2004 UTC (20 years ago) by gezelter
File size: 2127 byte(s)
Log Message:
Adding visualization parts

File Contents

# Content
1 #include <iostream>
2 #include <fstream>
3 #include <string>
4 #include <vector>
5 #include "visualizerCmd.h"
6
7
8 int count_tokens(char *line, char *delimiters);
9 using namespace std;
10
11 int main(int argc, char* argv[]){
12
13 gengetopt_args_info args_info;
14 int npts, i, j, k, write_size;
15 double xmin, xmax, x;
16 double ymin, ymax, y;
17 double zmin, zmax, z;
18 double r, theta, phi, f, s, w, w0;
19 double fp, sp, wp, tot;
20 double rup, ru, rlow;
21 FILE *my_file;
22
23 //parse the command line options
24 if (cmdline_parser (argc, argv, &args_info) != 0)
25 exit(1) ;
26
27 if (args_info.input_given){
28 shapeFileName = args_info.input_arg;
29 }
30 else{
31 std::cerr << "Does not have shape file name" << endl;
32 exit(1);
33 }
34
35 // grid has a default value (default=51), so it is always given
36 npts = args_info.grid_arg;
37
38 printf ("opening %s\n", shapeFileName);
39 shapeFile = fopen( shapeFileName, "r" );
40
41
42
43 my_file = fopen("junk.dat", "w");
44 if (my_file == NULL) {
45 (void) fprintf(stderr, "No File\n");
46 exit(8);
47 }
48
49 npts = 51;
50
51 xmin = -5.0;
52 xmax = 5.0;
53
54 ymin = -5.0;
55 ymax = 5.0;
56
57 zmin = -5.0;
58 zmax = 5.0;
59
60 rup = 4.0;
61 ru = 3.35;
62 rlow = 2.75;
63 w0 = 0.07715;
64
65 for (i = 0; i < npts; i++) {
66 x = xmin + (xmax-xmin) * (double)i/(double)npts;
67
68 for (j = 0; j < npts; j++) {
69 y = ymin + (ymax-ymin) * (double)j/(double)npts;
70
71 for (k = 0; k < npts; k++) {
72 z = zmin + (zmax-zmin) * (double)k/(double)npts;
73
74 r = sqrt(x*x + y*y + z*z);
75 theta = acos(z/r);
76 phi = atan(y/x);
77
78 fp = pow((cos(theta) - 0.6),2) * pow((cos(theta) + 0.8),2);
79
80 sp = (rup + 2.0*r - 3.0*rlow) * pow((rup - r),2) / pow((rup -rlow),3);
81
82 if (r > rup) sp = 0.0;
83 if (r < rlow) sp = 1.0;
84
85 wp = sp*(fp - w0);
86
87 f = sin(theta)*sin(2.0*theta)*cos(2.0*phi);
88
89 s = (ru + 2.0*r - 3.0*rlow) * pow((ru - r),2) / pow((ru -rlow),3);
90
91 if (r > ru) s = 0.0;
92 if (r < rlow) s = 1.0;
93
94 w = f*s;
95
96 tot = w + wp;
97
98 write_size = fwrite(&tot,1,8,my_file);
99
100 }
101 }
102 }
103 fclose(my_file);
104 return(0);
105 }