1 /*  XNECVIEW - a program for visualizing NEC2 input and output data
2  *
3  *  Copyright (C) 1998-2006, Pieter-Tjerk de Boer -- pa3fwm@amsat.org
4  *
5  *  Distributed on the conditions of version 2 of the GPL: see the files
6  *  README and COPYING, which accompany this source file.
7  *
8  *  This module contains almost all 3D drawing related stuff, like coordinate
9  *  transformation and drawing routines based on generic low-level graphics
10  *  routines, as provided by the Outdev-struct pointed to by *out.
11  *  The only 3D drawing code not included here, is the code for drawing the
12  *  the gain pattern as an opaque surface, i.e., with hidden-line removal;
13  *  that code is in draw_opaque.c, because of its size and complexity.
14  *
15  *  Note: the code for drawing the 2D plots of several quantities versus
16  *  frequency is in 'freqplot.c'.
17  *
18  */
19 
20 #include <stdio.h>
21 #include <math.h>
22 #include <stdlib.h>
23 #include <string.h>
24 #include <gdk/gdk.h>
25 
26 #include "xnecview.h"
27 
28 /* ------------------------- coordinate conversion and related arithmetic -------------------------------------------------- */
29 
30                                /* these define the projection: */
31 double phi=ini_phi,theta=ini_theta;    /* angles defining direction of eye */
32 double zoom=ini_zoom;                  /* zoom factor */
33 double trx=ini_trx,try=ini_try;        /* 2D-translation, as a fraction of winsize */
34 int winsizex,winsizey;                 /* size of window in pixels */
35 
36                                /* and these are calculated from them: */
37 double Xx,Xy,Xz,Yx,Yy,Yz;         /* projection matrix */
38 double Xtr,Ytr;                   /* 2D-translation */
39 int winsize;                      /* min(winsizex,winsizey) */
40 Point eyevec;                     /* tip of arrow pointing toward observer */
41 double Zx,Zy,Zz;                  /* weighting factors for calculating distance to observer (possibly in "locked" direction) */
42 
43 
44 
calcproj(void)45 void calcproj(void)     /* calculate the projection matrix etc. */
46 {
47    double ph,th;
48    double scale;
49 
50    if (winsizey<winsizex) winsize=winsizey;
51    else winsize=winsizex;
52    scale=winsize*zoom/(3.0*extr);
53    ph=phi*(M_PI/180.);
54    th=theta*(M_PI/180.);
55    Xz = 0;
56    Yz = -scale*sin(th);
57    Xx = -scale*sin(ph);
58    Xy = scale*cos(ph);
59    Yx = scale*cos(ph)*cos(th);
60    Yy = scale*sin(ph)*cos(th);
61    Xtr = 0.5 + trx*zoom*winsize + 0.5*winsizex;
62    Ytr = 0.5 + try*zoom*winsize + 0.5*winsizey;
63    eyevec.x = sin(th)*cos(ph);
64    eyevec.y = sin(th)*sin(ph);
65    eyevec.z = cos(th);
66    if (phasedirlock) {
67       ph=phasephi*(M_PI/180.);
68       th=phasetheta*(M_PI/180.);
69    }
70    Zz = cos(th);
71    Zx = sin(th)*cos(ph);
72    Zy = sin(th)*sin(ph);
73 }
74 
75 
proj(Point * p,double * c)76 void proj(Point *p,double *c)      /* perform the projection: *p defines a point in 3D-space, *c defines a point in the window */
77 {
78    c[0] = Xtr + Xx*p->x + Xy*p->y + Xz*p->z;
79    c[1] = Ytr + Yx*p->x + Yy*p->y + Yz*p->z;
80 }
81 
82 
83 
interpolate(double * p,double f,double * a,double * b)84 void interpolate(double *p,double f,double *a,double *b)   /* interpolate linearly between two points in the window (f=0..1 -> *p=*a..*b) */
85 {
86    p[0]=(1-f)*a[0]+f*b[0];
87    p[1]=(1-f)*a[1]+f*b[1];
88 }
89 
90 
91 
92 /* test whether a line through two points is "ascending" ( / ) as opposed to "descending" ( \ ) */
93 #define ascending(c1,c2) (((c2)[0]-(c1)[0])*((c2)[1]-(c1)[1])<0)
94 
95 /* ------------------------- antenna structure drawing ------------------------------------------------------------------- */
96 
draw_wires(void)97 void draw_wires(void)                /* draw the wires of the antenna */
98 {
99    Wire *wi;
100    double c0[2],c1[2];
101    char s[20];
102 
103    if (numwires==0) return;
104 
105    SetForeground(&c_wire);
106    SetLineAttributes(2, GDK_LINE_SOLID, GDK_CAP_ROUND, GDK_JOIN_ROUND);
107    for (wi=wires;wi<wires+numwires;wi++) {
108       proj(&wi->p0,c0);
109       proj(&wi->p1,c1);
110       DrawLine(c0[0], c0[1], c1[0], c1[1]);
111       if (structplot==SPtags && !quick && wi->itg>0 && wi->dispnr) {
112          sprintf(s,"%i",wi->itg);
113          if (ascending(c0,c1))
114             DrawStringUL((c0[0]+c1[0])/2+2,(c0[1]+c1[1])/2,s);
115          else
116             DrawStringLL((c0[0]+c1[0])/2+2,(c0[1]+c1[1])/2,s);
117       }
118    }
119 }
120 
121 
draw_loads(void)122 void draw_loads(void)                 /* (re)draw loaded wires */
123 {
124    Load *lo;
125    double c0[2],c1[2],c2[2],c3[2];
126 
127    if (numloads==0) return;
128    SetForeground(&c_load);
129    SetLineAttributes(3, GDK_LINE_SOLID, GDK_CAP_ROUND, GDK_JOIN_ROUND);
130    for (lo=loads;lo<loads+numloads;lo++) {
131       proj(&lo->wire->p0,c0);
132       proj(&lo->wire->p1,c1);
133       if (lo->wire->ns>0) {
134          interpolate(c2,(double)(lo->firstseg-1)/(lo->wire->ns),c0,c1);
135          interpolate(c3,(double)lo->lastseg/(lo->wire->ns),c0,c1);
136          DrawLine(c2[0], c2[1], c3[0], c3[1]);
137       } else {
138          DrawLine(c0[0], c0[1], c1[0], c1[1]);
139       }
140    }
141 }
142 
proj_segment(Wire * wire,int seg,double * c0,double * c1)143 void proj_segment(Wire *wire,int seg,double *c0,double *c1)
144 {
145    if (wire->ns>0) {
146       double d0[2],d1[2];
147       proj(&wire->p0,d0);
148       proj(&wire->p1,d1);
149       interpolate(c0,(double)(seg-1)/(wire->ns),d0,d1);
150       interpolate(c1,(double)seg/(wire->ns),d0,d1);
151    } else {
152       proj(&wire->p0,c0);
153       proj(&wire->p1,c1);
154    }
155 }
156 
157 
draw_excis(void)158 void draw_excis(void)                 /* (re)draw wires with excitations */
159 {
160    Exci *ex;
161    double c0[2],c1[2];
162 
163    if (numexcis==0) return;
164    SetForeground(&c_exci);
165    SetLineAttributes(4, GDK_LINE_SOLID, GDK_CAP_ROUND, GDK_JOIN_ROUND);
166    for (ex=excis;ex<excis+numexcis;ex++) {
167       proj_segment(ex->wire, ex->seg, c0, c1);
168       DrawLine(c0[0], c0[1], c1[0], c1[1]);
169    }
170 }
171 
172 
draw_netws(void)173 void draw_netws(void)                 /* (re)draw networks and transmission lines */
174 {
175    Netw *ne;
176    double c10[2],c11[2],c20[2],c21[2];  /* endpoints of relevant wires */
177 
178    if (numnetws==0) return;
179    SetForeground(&c_netw);
180    for (ne=netws;ne<netws+numnetws;ne++) {
181       proj_segment(ne->wire1, ne->seg1, c10, c11);
182       proj_segment(ne->wire2, ne->seg2, c20, c21);
183       if (ne->type==1) {
184          SetLineAttributes(1, GDK_LINE_SOLID, GDK_CAP_ROUND, GDK_JOIN_ROUND);
185          DrawLine(c10[0], c10[1], c20[0], c20[1]);
186          DrawLine(c11[0], c11[1], c21[0], c21[1]);
187       } else if (ne->type==-1) {
188          SetLineAttributes(1, GDK_LINE_SOLID, GDK_CAP_ROUND, GDK_JOIN_ROUND);
189          DrawLine(c10[0], c10[1], c21[0], c21[1]);
190          DrawLine(c11[0], c11[1], c20[0], c20[1]);
191       } else {
192          double d1[2],d2[2];
193          interpolate(d1,0.5,c10,c11);
194          interpolate(d2,0.5,c20,c21);
195          SetLineAttributes(2, GDK_LINE_SOLID, GDK_CAP_ROUND, GDK_JOIN_ROUND);
196          DrawLine(d1[0], d1[1], d2[0], d2[1]);
197       }
198    }
199 }
200 
201 
draw_surfaces(void)202 void draw_surfaces(void)
203 {
204    Surface *su;
205    double c0[2],c1[2],c2[2],c3[2];
206    double cc[2];
207 
208    if (numsurfaces==0) return;
209    SetLineAttributes(0, GDK_LINE_SOLID, GDK_CAP_ROUND, GDK_JOIN_ROUND);
210    for (su=surfaces;su<surfaces+numsurfaces;su++) {
211       if ( eyevec.x * (su->pa.x-su->pc.x) + eyevec.y * (su->pa.y-su->pc.y) + eyevec.z * (su->pa.z-su->pc.z) >= 0 )
212          SetForeground(&c_surf);
213       else
214          SetForeground(&c_back);
215       proj(&su->pc,cc);
216       proj(&su->p0,c0);
217       proj(&su->p1,c1);
218       switch (su->type) {
219          case SU_rect:
220          case SU_arb:
221             c2[0]=2*cc[0]-c0[0]; c2[1]=2*cc[1]-c0[1];
222             c3[0]=2*cc[0]-c1[0]; c3[1]=2*cc[1]-c1[1];
223             break;
224          case SU_quad:
225             proj(&su->p3,c3);
226             /* fallthrough */
227          case SU_tri:
228             proj(&su->p2,c2);
229             break;
230       }
231       if (su->type==SU_tri) {
232          /* triangle; draw border: */
233          DrawLine(c0[0], c0[1], c1[0], c1[1]);
234          DrawLine(c1[0], c1[1], c2[0], c2[1]);
235          DrawLine(c2[0], c2[1], c0[0], c0[1]);
236       } else {
237          /* joint drawing code for all other shapes, which are drawn as various quadrilaterals */
238          /* draw diagonals: */
239          DrawLine(c0[0], c0[1], c2[0], c2[1]);
240          DrawLine(c1[0], c1[1], c3[0], c3[1]);
241          if (su->type!=SU_arb) {
242             /* draw border: */
243             DrawLine(c0[0], c0[1], c1[0], c1[1]);
244             DrawLine(c1[0], c1[1], c2[0], c2[1]);
245             DrawLine(c2[0], c2[1], c3[0], c3[1]);
246             DrawLine(c3[0], c3[1], c0[0], c0[1]);
247          } else {
248             /* draw additional diagonals */
249             DrawLine(0.7071*(c0[0]+c1[0])-0.4142*cc[0], 0.7071*(c0[1]+c1[1])-0.4142*cc[1], 0.7071*(c2[0]+c3[0])-0.4142*cc[0], 0.7071*(c2[1]+c3[1])-0.4142*cc[1]);
250             DrawLine(0.7071*(c0[0]+c3[0])-0.4142*cc[0], 0.7071*(c0[1]+c3[1])-0.4142*cc[1], 0.7071*(c2[0]+c1[0])-0.4142*cc[0], 0.7071*(c2[1]+c1[1])-0.4142*cc[1]);
251          }
252       }
253    }
254 }
255 
256 
257 
draw_antenna(int ie)258 void draw_antenna(int ie)              /* draw the entire antenna */
259 {
260    draw_wires();
261    draw_loads();
262    draw_excis();
263    draw_netws();
264    if (Pending() && ie) return;
265    draw_surfaces();
266 }
267 
268 
269 /* ------------------------- currents/charges drawing -------------------------------------------------------------------- */
270 
cmpsegdist(const Segcur ** s1,const Segcur ** s2)271 int cmpsegdist(const Segcur **s1, const Segcur **s2)
272 {
273    double d;
274    d= (*s1)->dist - (*s2)->dist;
275    if (d<0) return -1;
276    if (d>0) return 1;
277    return 0;
278 }
279 
280 
addvector(Point * pe,float a[3],double q)281 void addvector(Point *pe,float a[3],double q)
282 {
283    pe->x+=q*a[0];
284    pe->y+=q*a[1];
285    pe->z+=q*a[2];
286 }
287 
draw_currents_anim(int ie,Currents * cu)288 void draw_currents_anim(int ie,Currents *cu)      /* draw the currents and charge distribution as vectors and circles with animation */
289 {
290    Segcur *se;
291    double si,co;
292    double sc;
293    double l;
294 
295    if (cu==NULL || cu->numseg==0) return;
296    l=0;
297    l+=(cu->s[0].p1.x-cu->s[0].p0.x)*(cu->s[0].p1.x-cu->s[0].p0.x);
298    l+=(cu->s[0].p1.y-cu->s[0].p0.y)*(cu->s[0].p1.y-cu->s[0].p0.y);
299    l+=(cu->s[0].p1.z-cu->s[0].p0.z)*(cu->s[0].p1.z-cu->s[0].p0.z);
300    l=sqrt(l);
301 
302    si = sin(-animphase);
303    co = cos(-animphase);
304    sc = iscale * 0.8 * l / cu->maxI;
305 
306    for (se=cu->s;se<cu->s+cu->numanimseg;se++) {
307       double c0[2],c1[2];
308       double q;
309       Point pe;
310       double ev[3];
311       int i;
312 
313       for (i=0;i<3;i++) ev[i] = se->re[i]*co + se->im[i]*si;
314       q = qscale * l * 0.5 * (se->qre*co + se->qim*si) / cu->maxQ;
315 
316       pe.x = se->c.x + sc*ev[0]; pe.y = se->c.y + sc*ev[1]; pe.z = se->c.z + sc*ev[2];
317       proj(&se->c,c0);
318       proj(&pe,c1);
319       SetLineAttributes(2, GDK_LINE_SOLID, GDK_CAP_ROUND, GDK_JOIN_ROUND);
320       SetForeground(&c_wire);
321       DrawLine(c0[0], c0[1], c1[0], c1[1]);
322 
323       SetLineAttributes(1, GDK_LINE_SOLID, GDK_CAP_ROUND, GDK_JOIN_ROUND);
324       if (q<0) SetForeground(&c_qneg); else SetForeground(&c_qpos);
325       q=fabs(q);
326       pe=se->c;
327       addvector(&pe,se->q0,q);  proj(&pe,c0);
328       addvector(&pe,se->q0,-q); addvector(&pe,se->q1,q);  proj(&pe,c1); DrawLine(c0[0], c0[1], c1[0], c1[1]);
329       addvector(&pe,se->q0,-q); addvector(&pe,se->q1,-q); proj(&pe,c0); DrawLine(c0[0], c0[1], c1[0], c1[1]);
330       addvector(&pe,se->q0,q);  addvector(&pe,se->q1,-q); proj(&pe,c1); DrawLine(c0[0], c0[1], c1[0], c1[1]);
331       addvector(&pe,se->q0,q);  addvector(&pe,se->q1,q);  proj(&pe,c0); DrawLine(c0[0], c0[1], c1[0], c1[1]);
332    }
333 }
334 
335 
draw_currents_noanim(int ie,Currents * cu)336 void draw_currents_noanim(int ie,Currents *cu)      /* draw the currents distribution as a static phase/magnitude display */
337 {
338    Segcur *se;
339    Segcur **sse;
340    int i;
341    double c0[2],c1[2];
342    double a=0,phase;
343    double dx,dy,l;
344    double w;
345 #ifdef GAINCALC
346    double totalr=0,totali=0;
347 #endif
348 
349    if (cu==NULL || cu->numseg==0) return;
350    sse=malloc(cu->numseg*sizeof(Segcur*));
351    if (sse==NULL) return;
352 
353    /* first, calculate the distance of each segment (defined in xnecview.h),
354       and at the same time fill the array **sse that will be sorted next */
355    se=cu->s;
356    for (i=0;i<cu->numseg;i++) {
357       se->dist=Zx*se->c.x+Zy*se->c.y+Zz*se->c.z;
358       sse[i]=se;
359       se++;
360    }
361 
362    /* next, sort the array elements according to this distance, so the element closest to the observer will be drawn last */
363    qsort(sse, cu->numseg, sizeof(Segcur*),
364       (int (*)(const void *, const void *))cmpsegdist);
365 
366    w = 360/(299.8/neco[rp_index].f);
367    for (i=0;i<cu->numseg;i++) {
368       se=sse[i];
369       proj(&se->p0,c0);
370       proj(&se->p1,c1);
371 
372       dx=c1[0]-c0[0];
373       dy=c1[1]-c0[1];
374       if (!quick) {
375          l=sqrt(dx*dx+dy*dy);
376          if (l==0) continue;
377          phase=se->phase;
378          switch (polarization) {
379             case POLnone:
380             case POLcolour:
381                a=se->a;
382                break;
383             case POLhor:
384                a=se->a*dx/l;
385                break;
386             case POLvert:
387                a=se->a*dy/l;
388                break;
389             case POLlhcp:
390                a=se->a;
391                phase+=180/M_PI*atan2(dy,dx);
392                break;
393             case POLrhcp:
394                a=se->a;
395                phase-=180/M_PI*atan2(dy,dx);
396                break;
397          }
398          if (a<0) { phase+=180; a=-a; }
399          if (distcor) phase+=w*se->dist;
400          phase+=phaseoffset;
401          while (phase<0) phase+=360;
402          while (phase>=360) phase-=360;
403       }
404 
405       if (!quick && maxthickness*a>=0.5) {
406          SetLineAttributes(maxthickness*a+0.5, GDK_LINE_SOLID, GDK_CAP_ROUND, GDK_JOIN_ROUND);
407          SetForeground(c_currents+(int)(phase*NC_phase/360));
408       } else {
409          SetLineAttributes(1, GDK_LINE_SOLID, GDK_CAP_ROUND, GDK_JOIN_ROUND);
410          SetForeground(&c_inactive);
411       }
412       DrawLine(c0[0], c0[1], c1[0], c1[1]);
413 
414 #ifdef GAINCALC
415       totalr+=a*l*cos(phase*M_PI/180);
416       totali+=a*l*sin(phase*M_PI/180);
417 #endif
418    }
419 #ifdef GAINCALC
420    totalr = sqrt(totalr*totalr+totali*totali);
421    printf("-->  %12g %12g\n",totalr,20*log10(totalr));
422 #endif
423 
424    free(sse);
425 }
426 
draw_phasecircle(void)427 void draw_phasecircle(void)
428 {
429    int xc,yc;  /* coordinates of center */
430    int r;      /* radius */
431    int i,j,x,y;
432    double phase,step;
433    double rc,rs;
434    double sintab[NC_phase/4];
435    double costab[NC_phase/4];
436 
437    r=winsize/10;
438    xc=yc=r+r/2+1;
439    SetLineAttributes(0, GDK_LINE_SOLID, GDK_CAP_ROUND, GDK_JOIN_ROUND);
440    step = 1/(2*M_PI*r) * 360;
441 
442    for (i=0;i<NC_phase/4;i++) {
443       phase=((i*360.0/NC_phase)-phaseoffset)*M_PI/180.0;
444       sintab[i]=sin(phase);
445       costab[i]=cos(phase);
446    }
447 
448    j=3*(int)(2*M_PI*r/NC_phase+2);
449    do {
450       j/=3;
451       SetLineAttributes(j, GDK_LINE_SOLID, GDK_CAP_ROUND, GDK_JOIN_ROUND);
452       for (i=0;i<NC_phase/4;i++) {
453          rs=r*sintab[i];
454          rc=r*costab[i];
455          x=xc+rc;
456          y=yc+rs;
457          SetForeground(c_currents+i);
458          DrawLine(xc,yc,x,y);
459          x=xc-rs;
460          y=yc+rc;
461          SetForeground(c_currents+i+NC_phase/4);
462          DrawLine(xc,yc,x,y);
463          x=xc-rc;
464          y=yc-rs;
465          SetForeground(c_currents+i+NC_phase/2);
466          DrawLine(xc,yc,x,y);
467          x=xc+rs;
468          y=yc-rc;
469          SetForeground(c_currents+i+3*NC_phase/4);
470          DrawLine(xc,yc,x,y);
471       }
472    } while (j>0);
473 }
474 
475 
476 /* ------------------------- near field drawing -------------------------------------------------------------------------- */
477 
draw_nearfield(int ie,NECoutput * ne)478 void draw_nearfield(int ie,NECoutput *ne)      /* draw the E and H field and Poynting vectors */
479 {
480    double f;
481    Nearfield *nf;
482    double sce,sch,scp;
483    double si,co;
484 
485    nf = ne->nf;
486    if (nf==NULL) return;
487 
488    if (nf->next==NULL) {
489       f = 1;
490    } else {
491       f = (nf->next->p.x-nf->p.x)*(nf->next->p.x-nf->p.x);
492       f += (nf->next->p.y-nf->p.y)*(nf->next->p.y-nf->p.y);
493       f += (nf->next->p.z-nf->p.z)*(nf->next->p.z-nf->p.z);
494       f = sqrt(f);
495    }
496    sce=f*escale/ne->maxe;
497    sch=f*hscale/ne->maxh;
498    scp=2*f*escale*hscale/(ne->maxe*ne->maxh);
499 
500    SetLineAttributes(0, GDK_LINE_SOLID, GDK_CAP_ROUND, GDK_JOIN_ROUND);
501    si = sin(-animphase);
502    co = cos(-animphase);
503    while (nf) {
504       Point pe;
505       double ev[3], hv[3], pv[3];
506       double c0[2],c1[2];
507       int valid;
508 
509       proj(&nf->p,c0);
510       valid=1;
511 
512       if (nf->evalid) {
513          int i;
514          for (i=0;i<3;i++) ev[i] = nf->ere[i]*co + nf->eim[i]*si;
515          pe.x = nf->p.x + sce*ev[0]; pe.y = nf->p.y + sce*ev[1]; pe.z = nf->p.z + sce*ev[2];
516          proj(&pe,c1);
517          SetForeground(&c_efield);
518          DrawLine(c0[0], c0[1], c1[0], c1[1]);
519       } else valid=0;
520 
521       if (nf->hvalid) {
522          int i;
523          for (i=0;i<3;i++) hv[i] = nf->hre[i]*co + nf->him[i]*si;
524          pe.x = nf->p.x + sch*hv[0]; pe.y = nf->p.y + sch*hv[1]; pe.z = nf->p.z + sch*hv[2];
525          proj(&pe,c1);
526          SetForeground(&c_hfield);
527          DrawLine(c0[0], c0[1], c1[0], c1[1]);
528       } else valid=0;
529 
530       if (valid && show_p) {
531          pv[0] = ev[1]*hv[2] - ev[2]*hv[1];
532          pv[1] = ev[2]*hv[0] - ev[0]*hv[2];
533          pv[2] = ev[0]*hv[1] - ev[1]*hv[0];
534          pe.x = nf->p.x + scp*pv[0]; pe.y = nf->p.y + scp*pv[1]; pe.z = nf->p.z + scp*pv[2];
535          proj(&pe,c1);
536          SetForeground(&c_poynting);
537          DrawLine(c0[0], c0[1], c1[0], c1[1]);
538       }
539 
540       nf=nf->next;
541    }
542 }
543 
544 
545 
546 /* ------------------------- gain pattern drawing (non-opaque) ----------------------------------------------------------- */
547 
548 
draw_polref(void)549 void draw_polref(void)
550 {
551    SetForeground(&c_gain_lin); DrawString(winsizex-5,5+0*fontheight,"linear",1,1);
552    SetForeground(&c_gain_lhcp); DrawString(winsizex-5,5+1*fontheight,"left-hand circular",1,1);
553    SetForeground(&c_gain_rhcp); DrawString(winsizex-5,5+2*fontheight,"right-hand circular",1,1);
554 }
555 
556 
setpolcolour(double axial)557 void setpolcolour(double axial)
558 {
559    if (axial>Polthr) SetForeground(&c_gain_lhcp);
560    else if (axial<-Polthr) SetForeground(&c_gain_rhcp);
561    else SetForeground(&c_gain_lin);
562 }
563 
564 
draw_gain_theta(int i,Radpattern * rp)565 void draw_gain_theta(int i,Radpattern *rp)
566 {
567    int j;
568    double c0[2],c1[2];
569    double a0,a1;
570 
571    proj(&rp->gpo[0][i],c1);
572    a1=rp->gain[0][i].axial;
573    for (j=1;j<rp->numphi;j++) {
574       c0[0]=c1[0]; c0[1]=c1[1];
575       proj(&rp->gpo[j][i],c1);
576       if (polarization==POLcolour) {
577          a0=a1;
578          a1=rp->gain[j][i].axial;
579          setpolcolour((a0+a1)/2);
580       }
581       DrawLine(c0[0], c0[1], c1[0], c1[1]);
582    }
583    c0[0]=c1[0]; c0[1]=c1[1];
584    proj(&rp->gpo[0][i],c1);
585    if (polarization==POLcolour) {
586       a0=a1;
587       a1=rp->gain[0][i].axial;
588       setpolcolour((a0+a1)/2);
589    }
590    DrawLine(c0[0], c0[1], c1[0], c1[1]);
591 }
592 
draw_gain_phi(int j,Radpattern * rp)593 void draw_gain_phi(int j,Radpattern *rp)
594 {
595    int i;
596    double c0[2],c1[2];
597    double a0,a1;
598 
599    proj(&rp->gpo[j][0],c1);
600    a1=rp->gain[j][0].axial;
601    for (i=1;i<rp->numtheta;i++) {
602       c0[0]=c1[0]; c0[1]=c1[1];
603       proj(&rp->gpo[j][i],c1);
604       if (polarization==POLcolour) {
605          a0=a1;
606          a1=rp->gain[j][i].axial;
607          setpolcolour((a0+a1)/2);
608       }
609       DrawLine(c0[0], c0[1], c1[0], c1[1]);
610    }
611 }
612 
613 
614 #define R90 (M_PI/2)
615 #define R180 M_PI
616 #define R270 (1.5*M_PI)
617 #define R360 (2*M_PI)
618 
619 #include <sys/time.h>
620 
draw_gain(int ie,Radpattern * rp)621 void draw_gain(int ie,Radpattern *rp)
622 {
623    int i,j;
624    int gainplot2=gainplot;
625 
626    SetForeground(&c_gain);
627    SetLineAttributes(0, GDK_LINE_SOLID, GDK_CAP_ROUND, GDK_JOIN_ROUND);
628 
629    if (gainplot==GPopaque && !quick) {
630       Radpattern *p;
631       int ok=0;
632       p=rp;
633       while (p) {
634          if (!opaque_impossible(p)) {
635             draw_opaque(p);
636             ok=1;
637          }
638          p=p->next;
639       }
640       if (ok) return;
641       fprintf(stderr,"Opaque gain plot not possible; reason: %s.\n",opaque_impossible(rp));
642       gainplot2=GPframe;
643    }
644 
645    while (rp) {
646       if (gainplot2==GPframe && !quick) {
647 
648          for (i=0;i<rp->numtheta;i++) draw_gain_theta(i,rp);
649          if (Pending() && ie) return;
650          for (j=0;j<rp->numphi;j++) draw_gain_phi(j,rp);
651 
652       } else {
653 
654          if (!scaleplot || theta!=90)      /* suppress if scaleplotting and XY-plane is perpendicular to screen */
655             for (i=0;i<rp->numtheta;i++)
656                if (fabs(rp->gtheta[i]-R90)<0.001)     /* points in XY-plane? */
657                   draw_gain_theta(i,rp);
658 
659          for (j=0;j<rp->numphi;j++)
660             if ((
661                   ( (fabs((rp->gphi[j]+R90)-R90)<0.001 || fabs((rp->gphi[j])-R180)<0.001)  /* points in XZ-plane? */
662                      && (!scaleplot || (phi!=0 && theta!=0)) )  /* suppress if that plane is perpendicular to screen and we're plotting the scale */
663                ) || (
664                   ( fabs((rp->gphi[j])-R90)<0.001 || fabs((rp->gphi[j])-R270)<0.001)       /* points in YZ-plane? */
665                      && (!scaleplot || (phi!=90 && theta!=0) )  /* suppress if that plane is perpendicular to screen and we're plotting the scale */
666                ))
667                draw_gain_phi(j,rp);
668 
669       }
670       rp=rp->next;
671    }
672 }
673 
674 
draw_gainscale(int ie)675 void draw_gainscale(int ie)
676 {
677    double a;
678    Point po;
679    double c0[2],c1[2],orig[2];
680 #define Npt 90
681 #define da R360/Npt
682    static double si[Npt],co[Npt];
683    static int firsttime=1;
684    static double radii[numGS][6]=     /* this array determines the radius of all gain scale circles; note: if you change this, also update the labels in GSnumbers[][] and in labels, here below */
685            {
686              {    1,     0.794328,    0.501187,    0.251189,    0.1,       0   },   /* lin.power:    0, -1, -3, -6, -10 dB */
687              {    1,     0.891251,    0.707946,    0.501187,    0.316228,  0   },   /* lin.voltage:  0, -1, -3, -6, -10 dB */
688              {    1,     0.839624,    0.558406,    0.311817,    0.174121,  0   },   /* arrl:         0, -3, -10, -20, -30 dB */
689              {    1,     0.75,        0.5,         0.25,        0   }               /* log:          0, -10, -20, -30 dB */
690            } ;
691    static char labels[numGS][6][4]=
692            {
693              { "0dB", "-1", "-3", "-6", "-10" },
694              { "0dB", "-1", "-3", "-6", "-10" },
695              { "0dB", "-3", "-10", "-20", "-30" },
696              { "0dB", "-10", "-20", "-30" }
697            } ;
698    int i,j;
699    double r;
700    double *rp;
701 
702    if (!scaleplot) return;
703 
704    SetForeground(&c_scale);
705    SetLineAttributes(0, GDK_LINE_SOLID, GDK_CAP_ROUND, GDK_JOIN_ROUND);
706 
707    if (firsttime) {
708       for (i=Npt-1, a=0; i>=0; i--, a+=da) {
709          si[i]=sin(a);
710          co[i]=cos(a);
711       }
712       firsttime=0;
713    }
714 
715    if (phi!=90 && theta!=0) {
716       /* draw radial lines (yz plane, 10 degrees spacing) */
717       po.x=0; po.y=0; po.z=0;
718       proj(&po,orig);
719       for (i=0;i<360;i=i+10) {
720          po.x=0;
721          po.y=extr*cos(i*M_PI/180);
722          po.z=extr*sin(i*M_PI/180);
723          proj(&po,c1);
724          DrawLine(orig[0],orig[1],c1[0],c1[1]);
725       }
726       /* draw circles */
727       rp=radii[gainscale];
728       j=0;
729       while ((r=*rp++*extr)) {
730          po.x=0; po.y=r; po.z=0;
731          proj(&po,c0);
732          for (i=0;i<Npt;i++) {
733             po.y=r*co[i];
734             po.z=r*si[i];
735             proj(&po,c1);
736             DrawLine(c0[0], c0[1], c1[0], c1[1]);
737             c0[0]=c1[0]; c0[1]=c1[1];
738          }
739          po.y=0; po.z=r; proj(&po,c1);
740          DrawStringLL(c1[0],c1[1],labels[gainscale][j++]);
741       }
742       if (Pending() && ie) return;
743    }
744 
745    if (theta!=90) {
746       /* draw radial lines (xy plane, 10 degrees spacing) */
747       po.x=0; po.y=0; po.z=0;
748       proj(&po,orig);
749       for (i=0;i<360;i=i+10) {
750          po.z=0;
751          po.x=extr*cos(i*M_PI/180);
752          po.y=extr*sin(i*M_PI/180);
753          proj(&po,c1);
754          DrawLine(orig[0],orig[1],c1[0],c1[1]);
755       }
756       /* draw circles */
757       rp=radii[gainscale];
758       j=0;
759       while ((r=*rp++*extr)) {
760          po.x=r; po.y=0; po.z=0;
761          proj(&po,c0);
762          for (i=0;i<Npt;i++) {
763             po.x=r*co[i];
764             po.y=r*si[i];
765             proj(&po,c1);
766             DrawLine(c0[0], c0[1], c1[0], c1[1]);
767             c0[0]=c1[0]; c0[1]=c1[1];
768          }
769          po.x=-r; po.y=0; proj(&po,c1);
770          DrawStringLL(c1[0],c1[1],labels[gainscale][j++]);
771       }
772       if (Pending() && ie) return;
773    }
774 
775    if (phi!=0 && theta!=0) {
776       /* draw radial lines (xz plane, 10 degrees spacing) */
777       po.x=0; po.y=0; po.z=0;
778       proj(&po,orig);
779       for (i=0;i<360;i=i+10) {
780          po.y=0;
781          po.x=extr*cos(i*M_PI/180);
782          po.z=extr*sin(i*M_PI/180);
783          proj(&po,c1);
784          DrawLine(orig[0],orig[1],c1[0],c1[1]);
785       }
786       /* draw circles */
787       rp=radii[gainscale];
788       j=0;
789       while ((r=*rp++*extr)) {
790          po.x=r; po.y=0; po.z=0;
791          proj(&po,c0);
792          for (i=0;i<Npt;i++) {
793             po.x=r*co[i];
794             po.z=r*si[i];
795             proj(&po,c1);
796             DrawLine(c0[0], c0[1], c1[0], c1[1]);
797             c0[0]=c1[0]; c0[1]=c1[1];
798          }
799          po.x=0; po.z=r; proj(&po,c1);
800          DrawStringLL(c1[0],c1[1],labels[gainscale][j++]);
801       }
802    }
803 }
804 
805 
806 /* ------------------------- draw complete picture ------------------------------------------------------------------------- */
807 
draw_axes()808 void draw_axes()		/* draw the axes */
809 {
810    double c0[2],c1[2],c2[2],c3[2];
811    Point pp;
812 
813    SetForeground(&c_axis);
814    SetLineAttributes(0, GDK_LINE_SOLID, GDK_CAP_ROUND, GDK_JOIN_ROUND);
815    pp.x=0; pp.y=0; pp.z=0; proj(&pp,c0);
816    pp.x=Axislen*extr; pp.y=0; pp.z=0; proj(&pp,c1);
817    pp.x=0; pp.y=Axislen*extr; pp.z=0; proj(&pp,c2);
818    pp.x=0; pp.y=0; pp.z=Axislen*extr; proj(&pp,c3);
819    DrawLine(c0[0], c0[1], c1[0], c1[1]);
820    DrawLine(c0[0], c0[1], c2[0], c2[1]);
821    DrawLine(c0[0], c0[1], c3[0], c3[1]);
822    if (!quick) {
823       if (c1[0]>=c0[0] || !ascending(c0,c1)) DrawStringLL(c1[0]+2,c1[1],"X"); else DrawStringUL(c1[0]+2,c1[1],"X");
824       if (c2[0]>=c0[0] || !ascending(c0,c2)) DrawStringLL(c2[0]+2,c2[1],"Y"); else DrawStringUL(c2[0]+2,c2[1],"Y");
825       if (c3[0]>=c0[0] || !ascending(c0,c3)) DrawStringLL(c3[0]+2,c3[1],"Z"); else DrawStringUL(c3[0]+2,c3[1],"Z");
826    }
827 
828    SetLineAttributes(0, GDK_LINE_ON_OFF_DASH, GDK_CAP_ROUND, GDK_JOIN_ROUND);
829    pp.x=0; pp.y=0; pp.z=0; proj(&pp,c0);
830    pp.x=-Axislen*extr; pp.y=0; pp.z=0; proj(&pp,c1);
831    DrawLine(c0[0], c0[1], c1[0], c1[1]);
832    pp.x=0; pp.y=-Axislen*extr; pp.z=0; proj(&pp,c1);
833    DrawLine(c0[0], c0[1], c1[0], c1[1]);
834    pp.x=0; pp.y=0; pp.z=-Axislen*extr; proj(&pp,c1);
835    DrawLine(c0[0], c0[1], c1[0], c1[1]);
836 }
837 
838 
draw_all(int ie)839 void draw_all(int ie)
840 {
841    char s[100];
842    do {
843       if (Pending() && ie) return;
844       ClearWindow();
845       if (gainplot==GPslice || gainplot==GPframe || gainplot==GPopaque) draw_gainscale(ie);
846       draw_axes();
847       if (Pending() && ie) break;
848       if (structplot) {
849          if (structplot==SPcurrents) {
850             if (rp_index>=0) {
851                draw_currents_noanim(ie,neco[rp_index].cu);
852                if (!quick) draw_phasecircle();
853             }
854          } else if (structplot==SPanim) {
855             if (rp_index>=0) {
856                draw_currents_anim(ie,neco[rp_index].cu);
857             }
858          }
859          else draw_antenna(ie);
860       }
861       if (Pending() && ie) break;
862       if (rp_index>=0) {
863          if (gainplot==GPslice || gainplot==GPframe || gainplot==GPopaque) {
864             draw_gain(ie,neco[rp_index].rp);
865             if (polarization==POLcolour) draw_polref();
866          }
867          if (gainplot==GPnearfield) draw_nearfield(ie,neco+rp_index);
868          if (neco[rp_index].rp) {
869             double maxgain,vgain;
870             if (polarization==POLnone || polarization==POLcolour) {
871                maxgain=neco[rp_index].d[neco_maxgain],
872                vgain=neco[rp_index].d[neco_vgain];
873             } else {
874                maxgain=neco[rp_index].d[Neco_polgain+Neco_gsize*polarization];
875                vgain=neco[rp_index].d[neco_vgain2];
876             }
877             sprintf(s,"f = %g MHz   maxgain = %g dBi   vgain = %g dBi",
878                neco[rp_index].f, maxgain, vgain);
879          } else {
880             sprintf(s,"f = %g MHz", neco[rp_index].f);
881 
882          }
883          SetForeground(&c_axis);
884          DrawString(winsizex/2,winsizey,s,0.5,0);
885       }
886    } while(0);
887    out->Complete();
888 }
889 
890 
891