@@ -1,282 +1,282 | |||
|
1 | 1 | #include <stdio.h> |
|
2 | 2 | #include <math.h> |
|
3 | 3 | #include <string.h> |
|
4 | 4 | #include <stdlib.h> |
|
5 | 5 | #include <fcntl.h> |
|
6 | 6 | #include <time.h> |
|
7 | 7 | #include "complex.h" |
|
8 | 8 | |
|
9 | 9 | #define NFREQ 512 |
|
10 | 10 | #define LNES 12 |
|
11 | 11 | #define LTES 16 |
|
12 | 12 | #define NANG 19 |
|
13 | 13 | |
|
14 | 14 | #define MAX(x,y) ((x)>(y) ? (x) : (y)) |
|
15 | 15 | #define MIN(x,y) ((x)<(y) ? (x) : (y)) |
|
16 | 16 | |
|
17 | 17 | fcomplex ****exlib; |
|
18 | 18 | int first=1; |
|
19 | 19 | |
|
20 | 20 | void swab_test(void); |
|
21 | 21 | //float swab(signed char *a); |
|
22 | 22 | void read_exlib(char *s,int nfreq, int lnes, int lTes, int nang, |
|
23 | 23 | fcomplex ****exlib); |
|
24 | 24 | |
|
25 | 25 | void initialize(void); |
|
26 | 26 | //void collision_(float *np, float *tp, float *fp, float *ap, float *yr, float *yi)// step 1 |
|
27 | 27 | void collision_(float *n, float *t, float *f, float *a, fcomplex *y); |
|
28 | 28 | |
|
29 | 29 | /* |
|
30 | 30 | void main(void){ |
|
31 | 31 | float n=4.9e11,t=1100.0,f=300.0,a=1.0; |
|
32 | 32 | int i; |
|
33 | 33 | fcomplex j; |
|
34 | 34 | for(i=0;i<10;i++){ |
|
35 | 35 | a+=0.05; |
|
36 | 36 | collision_(&n,&t,&f,&a,&j); |
|
37 | 37 | printf("%f %f %f \n",a,j.r,j.i); |
|
38 | 38 | } |
|
39 | 39 | } |
|
40 | 40 | */ |
|
41 | 41 | |
|
42 | 42 | //void collision_(float *np, float *tp, float *fp, float *ap, float *yr, float *yi){ // step 2 |
|
43 | 43 | void collision_(float *np, float *tp, float *fp, float *ap, fcomplex *y){ |
|
44 | 44 | float ne[]={5.0e10, 1.0e11, 2.0e11, 3.5e11, 5.0e11, 7.5e11, |
|
45 | 45 | 1.0e12, 1.5e12, 2.0e12, 2.5e12, 3.0e12, 3.5e12}; |
|
46 | 46 | float te[]={600.0, 700.0, 800.0, 900.0, 1000.0, 1150.0, 1300.0, 1500.0, |
|
47 | 47 | 1750.0, 2000.0, 2250.0, 2500.0, 2750.0, 3000.0, 3500.0, 4000.0}; |
|
48 | 48 | float alpha[]={.125, .1875, .25, .3125, .375, .5, .625, .75, 1.0, 1.25, |
|
49 | 49 | 1.5, 2.0, 2.5, 3.0, 3.5, 4.0, 4.5, 5.0, 6.0}; |
|
50 | 50 | float freq[2],fmax=6237.79; |
|
51 | 51 | float n,t,f,a; |
|
52 | 52 | float v,w1,w2,w3,w4,weight; |
|
53 | 53 | fcomplex g; |
|
54 | 54 | int i,j,k,l,m,i1,j1,k1,l1; |
|
55 | 55 | |
|
56 | 56 | if(first==1){ |
|
57 | 57 | //printf("initi"); |
|
58 | 58 | initialize(); |
|
59 | 59 | first=0; |
|
60 | 60 | } |
|
61 | 61 | //printf(" inputs in collision_ \n y.r %f y.i %f \n *************\n",(*y).r,(*y).i); |
|
62 | 62 | |
|
63 | 63 | //printf("After:print\n"); |
|
64 | 64 | a=*ap; t=*tp; f=*fp; n=*np; |
|
65 | 65 | |
|
66 | 66 | a=MAX(MIN(a,6.0),0.125); |
|
67 | 67 | t=MAX(MIN(t,4000.0),600.0); |
|
68 | 68 | f=MIN(f,fmax-1.0); |
|
69 | 69 | n=MAX(MIN(n,3.5e12),5.0e10); |
|
70 | 70 | |
|
71 | 71 | i1=((float)(NFREQ-1))*f/fmax; |
|
72 | 72 | i1=MIN(i1,NFREQ-1); |
|
73 | 73 | freq[0]=i1*fmax/(float)(NFREQ-1); |
|
74 | 74 | freq[1]=(i1+1)*fmax/(float)(NFREQ-1); |
|
75 | 75 | |
|
76 | 76 | j1=0; |
|
77 | 77 | while(n>ne[j1+1]) j1++; |
|
78 | 78 | j1=MIN(j1,LNES-1); |
|
79 | 79 | |
|
80 | 80 | k1=0; |
|
81 | 81 | while(t>te[k1+1]) k1++; |
|
82 | 82 | k1=MIN(k1,LTES); |
|
83 | 83 | |
|
84 | 84 | l1=0; |
|
85 | 85 | while(a>alpha[l1+1]) l1++; |
|
86 | 86 | l1=MIN(l1,NANG); |
|
87 | 87 | |
|
88 | 88 | // printf("i,j,k,l: %i %i %i %i \n",i1,j1,k1,l1); |
|
89 | 89 | |
|
90 | 90 | v=(fmax/(float)(NFREQ-1))* |
|
91 | 91 | (ne[j1+1]-ne[j1])*(te[k1+1]-te[k1])*(alpha[l1+1]-alpha[l1]); |
|
92 | 92 | g=Complex(0.0,0.0); |
|
93 | 93 | |
|
94 | 94 | for(i=0;i<2;i++){ |
|
95 | 95 | if(i==0) |
|
96 | 96 | w1=freq[1]-f; |
|
97 | 97 | else |
|
98 | 98 | w1=f-freq[0]; |
|
99 | 99 | for(j=0;j<2;j++){ |
|
100 | 100 | if(j==0) |
|
101 | 101 | w2=ne[j1+1]-n; |
|
102 | 102 | else |
|
103 | 103 | w2=n-ne[j1]; |
|
104 | 104 | for(k=0;k<2;k++){ |
|
105 | 105 | if(k==0) |
|
106 | 106 | w3=te[k1+1]-t; |
|
107 | 107 | else |
|
108 | 108 | w3=t-te[k1]; |
|
109 | 109 | for(l=0;l<2;l++){ |
|
110 | 110 | if(l==0) |
|
111 | 111 | w4=alpha[l1+1]-a; |
|
112 | 112 | else |
|
113 | 113 | w4=a-alpha[l1]; |
|
114 | 114 | weight=w1*w2*w3*w4; |
|
115 | 115 | g.r+=exlib[i1+i][j1+j][k1+k][l1+l].r*weight; |
|
116 | 116 | g.i+=exlib[i1+i][j1+j][k1+k][l1+l].i*weight; |
|
117 | 117 | // printf(" in collision i1+i %d j1+j %d k1+k %d l1+l %d exlib.r %f exlib.i %f \n",i1+i,j1+j,k1+k,l1+l,exlib[i1+i][j1+j][k1+k][l1+l].r,exlib[i1+i][j1+j][k1+k][l1+l].i); |
|
118 | 118 | // printf(" in collision g.r %f g.i %f \n",g.r,g.i); |
|
119 | 119 | } |
|
120 | 120 | } |
|
121 | 121 | } |
|
122 | 122 | } |
|
123 | 123 | g.r/=v; g.i/=v; |
|
124 | 124 | |
|
125 | 125 | *y=g; |
|
126 | 126 | //*yr=g.r |
|
127 | 127 | //*yi=g.i |
|
128 | 128 | //printf(" outputs in collision_ \n g.r %f g.i %f \n *************\n",g.r,g.i); |
|
129 | 129 | //exit(-2); |
|
130 | 130 | //getchar(); |
|
131 | 131 | return; |
|
132 | 132 | } |
|
133 | 133 | |
|
134 | 134 | void initialize(void){ |
|
135 | 135 | // char s[]="/usr/local/lib/faraday/jlib26feb2002"; |
|
136 | 136 | // char s[]="/usr/local/lib/faraday/jlib26feb2002"; |
|
137 | 137 | // FILE *s=fopen("/usr/local/lib/faraday/jlib26feb2002","r"); |
|
138 | 138 | int nfreq=NFREQ, lnes=LNES, lTes=LTES, nang=NANG; |
|
139 | 139 | int i,j,k; |
|
140 | 140 | int crr,crr2; |
|
141 | 141 | int the_number=1025; |
|
142 | 142 | |
|
143 | 143 | char the_path[1025]; |
|
144 | 144 | |
|
145 | 145 | get_path_reader_(the_path,&the_number); |
|
146 | 146 | //printf("C the_number: %i",the_number); |
|
147 | 147 | //printf("After\n"); |
|
148 | 148 | //printf("The Path: %s", the_path); |
|
149 | 149 | //printf("END\n"); |
|
150 | 150 | //printf("%s%n",the_path,&crr); |
|
151 | 151 | //printf("=%d",crr); |
|
152 | 152 | //printf("END2\n"); |
|
153 | 153 | char the_path_true[the_number]; |
|
154 | 154 | memcpy(the_path_true,the_path,the_number); |
|
155 | 155 | the_path_true[the_number] = '\0'; |
|
156 | printf("%s%n",the_path_true,&crr2); | |
|
156 | //printf("%s%n",the_path_true,&crr2); | |
|
157 | 157 | //printf("=%d",crr2); |
|
158 | 158 | //printf("END3\n"); |
|
159 | 159 | //printf("PATH: %s",s); |
|
160 | 160 | //fprintf(s); |
|
161 | 161 | exlib=malloc(nfreq*8); |
|
162 | 162 | for(i=0;i<nfreq;i++){ |
|
163 | 163 | exlib[i]=malloc(lnes*8); |
|
164 | 164 | for(j=0;j<lnes;j++){ |
|
165 | 165 | exlib[i][j]=malloc(lTes*8); |
|
166 | 166 | for(k=0;k<lTes;k++){ |
|
167 | 167 | exlib[i][j][k]=(fcomplex *)malloc(nang*sizeof(fcomplex)); |
|
168 | 168 | } |
|
169 | 169 | } |
|
170 | 170 | } |
|
171 | 171 | |
|
172 | 172 | read_exlib(the_path_true,nfreq,lnes,lTes,nang,exlib); |
|
173 | 173 | } |
|
174 | 174 | |
|
175 | 175 | float swab(signed char *a) |
|
176 | 176 | // table is in big endian format |
|
177 | 177 | { |
|
178 | 178 | float mant; |
|
179 | 179 | unsigned char *amant; |
|
180 | 180 | printf(" in swab 0 *** a[0] %c a[1] %c a[2] %c a[3] %c\n",a[0],a[1],a[2],a[3]); |
|
181 | 181 | amant=(unsigned char *) &mant; |
|
182 | 182 | amant[3]=a[0]; amant[2]=a[1]; amant[1]=a[2]; amant[0]=a[3]; |
|
183 | 183 | //printf(" in swab 1 ***\n"); |
|
184 | 184 | return(mant); |
|
185 | 185 | } |
|
186 | 186 | |
|
187 | 187 | void swab_test(void) |
|
188 | 188 | // table is in big endian format |
|
189 | 189 | { |
|
190 | 190 | printf("*********in swab stest ********\n"); |
|
191 | 191 | getchar(); |
|
192 | 192 | } |
|
193 | 193 | |
|
194 | 194 | |
|
195 | 195 | void read_exlib(char *the_path_true,int nfreq, int lnes, int lTes, int nang, |
|
196 | 196 | fcomplex ****exlib) |
|
197 | 197 | { |
|
198 | 198 | FILE *fp; |
|
199 | 199 | int i,j,k,l; |
|
200 | 200 | float a; |
|
201 | 201 | float t1, t2; |
|
202 | 202 | signed char *t1a,*t2b; |
|
203 | 203 | float mant; |
|
204 | 204 | unsigned char *amant; |
|
205 | 205 | printf("Start reading exlib file \n"); |
|
206 | 206 | if( (fp=fopen(the_path_true,"r")) == NULL){ |
|
207 | 207 | printf("Je Library file opening error\n"); |
|
208 | 208 | //printf("%s\n", strerror(errno)); |
|
209 | 209 | exit(-2); |
|
210 | 210 | } |
|
211 | 211 | for(l=0;l<nang;l++) |
|
212 | 212 | for(k=0;k<lTes;k++) |
|
213 | 213 | for(j=0;j<lnes;j++) |
|
214 | 214 | for(i=0;i<nfreq;i++){ |
|
215 | 215 | //printf(" in read_exlib 0***\n"); |
|
216 | 216 | fread(&t1,sizeof(float),1,fp); |
|
217 | 217 | fread(&t2,sizeof(float),1,fp); |
|
218 | 218 | //printf(" in read_exlib 1***\n"); |
|
219 | 219 | //printf(" i %d j %d k %d l %d t1 %f t2 %f \n",i,j,k,l,t1,t2); |
|
220 | 220 | t1a=(signed char *)&t1; |
|
221 | 221 | t2b=(signed char *)&t2; |
|
222 | 222 | //printf(" in read_exlib 1.5 t1a %d t2b %d ***\n",t1a,t2b); |
|
223 | 223 | ////t1=swab(t1a); |
|
224 | 224 | //printf(" t1 -- in readexlib 1.8 0 *** t1a[0] %c t1a[1] %c t1a[2] %c t1a[3] %c\n",t1a[0],t1a[1],t1a[2],t1a[3]); |
|
225 | 225 | amant=(unsigned char *) &mant; |
|
226 | 226 | amant[3]=t1a[0]; amant[2]=t1a[1]; amant[1]=t1a[2]; amant[0]=t1a[3]; |
|
227 | 227 | t1=mant; |
|
228 | 228 | |
|
229 | 229 | ////t2=swab(t2b); |
|
230 | 230 | //printf(" t2 -- in readexlib 1.8 0 *** t2b[0] %c t2b[1] %c t2b[2] %c t2b[3] %c\n",t2b[0],t2b[1],t2b[2],t2b[3]); |
|
231 | 231 | amant=(unsigned char *) &mant; |
|
232 | 232 | amant[3]=t2b[0]; amant[2]=t2b[1]; amant[1]=t2b[2]; amant[0]=t2b[3]; |
|
233 | 233 | t2=mant; |
|
234 | 234 | |
|
235 | 235 | ////t1=swab((signed char *)&t1); |
|
236 | 236 | ////t2=swab((signed char *)&t2); |
|
237 | 237 | //printf(" in read_exlib 2 i %d j %d k %d l %d t1 %f t2 %f \n",i,j,k,l,t1,t2); |
|
238 | 238 | exlib[i][j][k][l]=Complex(t1,t2); |
|
239 | 239 | //printf(" in read_exlib 3 i %d j %d k %d l %d t1 %f t2 %f \n",i,j,k,l,exlib[i][j][k][l].r,exlib[i][j][k][l].i); |
|
240 | 240 | //getchar(); |
|
241 | 241 | } |
|
242 | 242 | // getchar(); |
|
243 | 243 | |
|
244 | 244 | fclose(fp); |
|
245 | 245 | printf("Done reading exlib file \n"); |
|
246 | 246 | } |
|
247 | 247 | /* |
|
248 | 248 | void read_exlib(char *s,int nfreq, int lnes, int lTes, int nang, |
|
249 | 249 | fcomplex ****exlib) |
|
250 | 250 | { |
|
251 | 251 | FILE *fp; |
|
252 | 252 | int i,j,k,l; |
|
253 | 253 | float a; |
|
254 | 254 | float t1, t2; |
|
255 | 255 | |
|
256 | 256 | printf("Start reading exlib fileee \n"); |
|
257 | 257 | |
|
258 | 258 | if( (fp=fopen(s,"r")) == NULL){ |
|
259 | 259 | printf("Je Library file opening error\n"); |
|
260 | 260 | exit(-2); |
|
261 | 261 | } |
|
262 | 262 | |
|
263 | 263 | for(l=0;l<nang;l++) |
|
264 | 264 | for(k=0;k<lTes;k++) |
|
265 | 265 | for(j=0;j<lnes;j++) |
|
266 | 266 | for(i=0;i<nfreq;i++){ |
|
267 | 267 | printf(" in read_exlib 0***\n"); |
|
268 | 268 | fread(&t1,sizeof(float),1,fp); |
|
269 | 269 | fread(&t2,sizeof(float),1,fp); |
|
270 | 270 | printf(" in read_exlib 1***\n"); |
|
271 | 271 | printf(" i %d j %d k %d l %d t1 %f t2 %f \n",i,j,k,l,t1,t2); |
|
272 | 272 | t1=swab((char *)&t1); t2=swab((char *)&t2); |
|
273 | 273 | printf(" in read_exlib 2 i %d j %d k %d l %d t1 %f t2 %f \n",i,j,k,l,t1,t2); |
|
274 | 274 | exlib[i][j][k][l]=Complex(t1,t2); |
|
275 | 275 | printf(" in read_exlib 3 i %d j %d k %d l %d t1 %f t2 %f \n",i,j,k,l,exlib[i][j][k][l].r,exlib[i][j][k][l].i); |
|
276 | 276 | getchar(); |
|
277 | 277 | } |
|
278 | 278 | |
|
279 | 279 | fclose(fp); |
|
280 | 280 | printf("Done reading exlib file \n"); |
|
281 | 281 | } |
|
282 | 282 | */ |
General Comments 0
You need to be logged in to leave comments.
Login now