1 #ifdef HAVE_STDLIB_H
2 #include <stdlib.h>
3 #endif
4 #include <stdio.h>
5 #include <math.h>
6 #include "yagi.h"
7
8 /* This routine finds the self impedance of a dipole of arbitrary
9 length 'length' , radius 'r' at a wavelength 'wavelength'
10
11 Both the real part R and the imaginary part X are returned. R is
12 found from Antenna Theory, Analysis and Design, by Balans (Published
13 by Harper and Row).
14
15 The real part, Rin is found by computing the radiation resistance Rr
16 (equ 4-70, pp 124) then computing Rin=Rr/(sin(k l/2)^2, as given by
17 equ 7-30a, pp 294.
18
19 The reactive part, Xin is found by computing Xm (equ 7-33, pp 294)
20 then using xin=Xm/(sin(k l/2)^2 (equ 70-30b pp 294). */
21
self(double r,double length,double lambda,double * Rin,double * Xin)22 void self(double r, double length, double lambda, double *Rin, double *Xin)
23 {
24 double beta, mu, epsilon,eta,sin_bl,cos_bl,current_scale_factor;
25 double bl, ci_bl, si_bl, ci_2bl, si_2bl,Xm,Rr;
26
27 beta=2*M_PI/lambda;
28 mu=4*M_PI*1e-7;
29 epsilon=8.854187818e-12;
30 eta=sqrt(mu/epsilon);
31 bl=beta*length;
32
33 /* We often need si(bl) , si(2*bl), ci(bl), ci(2*bl). These are
34 time consuming to calculate, so I calculate just once. The Numerical
35 recipes routine cisi, calculates both si and ci at the same time,
36 so I wont use my routines ci and si, which both simply call cisi,
37 but was one of the results */
38 cisi(bl, &ci_bl, &si_bl);
39 cisi(2*bl, &ci_2bl, &si_2bl);
40 sin_bl=sin(bl);
41 cos_bl=cos(bl);
42 /* imaginary part of zin given by Balans, pp 294 */
43 Xm=(eta/(4*M_PI)) *
44 (2*si_bl+ cos_bl*( 2*si_bl - si_2bl )
45 -sin_bl*( 2*ci_bl-ci_2bl-ci(2*beta*r*r/length)) );
46 current_scale_factor=sin(bl/2)*sin(bl/2); /* sin(bl/2)^2 */
47 *Xin=Xm/current_scale_factor; /* Xin=Xm/(sin(bl/2)*sin(bl/2)); */
48 /* The radiation resistance, as described by Balans C. A. Antenna theory,
49 pp124, Published by Harper and Row, (1982) */
50
51 /* Rr=(eta/(M_PI*2.0))*(EULER+log(bl)-ci_bl+0.5*sin_bl*
52 (si_2bl-2.0*si_bl)+0.5*cos_bl
53 *(EULER+log(bl/2)+ci_2bl-2*ci_bl)); */
54 Rr=(eta/(2.0*M_PI))*(EULER+log(bl)-ci_bl+0.5*sin_bl*(si_2bl-2.0*si_bl)+0.5*cos_bl*(EULER+log(bl/2.0)+ci_2bl-2*ci_bl));
55 *Rin=Rr/current_scale_factor;
56 /* printf("Rr=%f Rin=%f csf=%f\n",Rr,*Rin,current_scale_factor); */
57 /* I'm a little bothered that Rin falls to ~60 Ohms at
58 resonance. Lawson suggested using 73, which would imply a large error.
59 A graph in RadCom by G3SEK would suggest the 60 Ohms is too low. A
60 graph in Krauss, would suggest ~70 Ohms, so what is going on????
61 I thought I'd use an equation in Krauss, but since this is the
62 Brown and King formuala given in Balans, that too gives ~60 Ohms.
63 The formula in Krass ignores the diameter, but gives much the same values
64 as the Balan's method */
65 /* Now compute by method in Krauss, pp 421. */
66 /* cot_bl_over_2=1.0/(tan(bl/2.0)); */
67
68 /* *Rin=30.0*( (1.0-(cot_bl_over_2*cot_bl_over_2))*Cin(2.0*bl)+4*cot_bl_over_2*cot_bl_over_2*Cin(bl)+2.0*cot_bl_over_2*(si_2bl-2.0*si_bl) );
69 printf("Rin=%f\n",*Rin); */
70 }
71