Logo Search packages:      
Sourcecode: yagiuda version File versions  Download package

yagi.c

#ifdef HAVE_STDLIB_H
#include <stdlib.h>
#endif

/* According to the method given in the book:
Yagi Antenna Design, by Dr. James L. Lawson, W2PV, ARRL.
Essentially, the self impedance Zxx of an elemnt is found. Then the mutural 
impedance between two elements Zxy is found. Then the currents and voltages
are computed  since:

I1 Z11 + I2 Z12 + I3 Z13 = V1
I1 Z21 + I2 Z22 + I3 Z23 = V2
I1 Z31 + I2 Z32 + I3 Z33 = V3

For just element 1 being driven, V2 = V3 = 0.
The input impedance is Zin= V_driven / I_driven
gain and other parameters can be found, as required. See book for
formulae */


#include <stdio.h>
#include <malloc.h>
#include <math.h>
#include <errno.h>
#include "yagi.h"

extern int optind, opterr, errno;

int main(int argc, char **argv)
{
      char *output_filename, *input_filename, *line;
      double **driven_data, **parasitic_data, **z, frequency, *v, **A, *b, *x;
      double d=0.0, real=0.0, imaginary=0.0, min_frequency=0.0, max_frequency=0.0;
      double step_frequency=0.0;  
      double angular_step=0.0;
      int   helpflg=0, errflg=0, pflg=0, dflg=0;
      int elements=-1, driven=-1, parasitic=-1, i=0, j=0, *indx;
      int element_number=0, c=0, sflg=0;
      FILE *ofp;
      struct FCOMPLEX *current, *voltage;
       /* necessary with gcc under DOS */
      opterr=0;

   while ((c =  getoptions(argc, argv, "dhps")) != -1)
   switch       (c) 
   {
                   case 'd': /* display currents */
                        dflg=1;
                        break;
                   case 'h':   /* help */
                        helpflg=1;
                        break;
                  case 'p':
                        pflg=1;
                        break;
                   case 's':   /* suppress all messages*/
                        sflg=1;
                        break;
                   case '?':
                        errflg=1;
                        break;
      }

      if(errflg || helpflg || (argc-optind)!=1)
      {
            usage_yagi(argv[0]);
            exit(1); 
      }
      /* allocate memory, using routines similar to that in Numerical Recipes */
      input_filename=string(0L,90L);
      line=string(0L,MAX_LINE);
      output_filename=string(0L,90L);
      output_filename=get_data_filenames(optind, argv, input_filename);
      if(errno)
         fprintf(stderr, "Errno = %d in yagi.c\n", errno);
      /* since the size of all the arrays depends on the number of elements,
      lets first find out how many elements we have,including driven and 
      parasitic sepparately.  Then try to allocate all the memory we need. If we 
      fail, then not much time is wasted. */
      elements=get_number_of_elements(input_filename, &driven, &parasitic);
      indx=ivector(1,2L*elements);

      /* allocate all memory */
      driven_data=dmatrix(1L,(long)(driven),1L,6L); /* 6 bits of info on each driven ele */
      if(parasitic > 0)
            parasitic_data=dmatrix(1L,(long) (parasitic),1L,4L); /* 4 on each parasitic ele */
      /* the impedance matrix for the impedance Z11 to Znn, needs an NxN
      matrix. Since the data is complex, it needs 2*N by N, where N is
      the total number of elelemnts */
      z=dmatrix(1L,(long) elements, 1L, elements*2L);
      A=dmatrix(1L, 2L*elements, 1L, 2L*elements);   /* A.x=b for large matrices */
      x=dvector(1L, 2L*elements);
      b=dvector(1L, 2L*elements);
      /* read text file containing yagi description into memory */
      read_yagi_data(line, input_filename, &frequency, &min_frequency, &max_frequency, 
      &step_frequency, driven, driven_data, parasitic, parasitic_data, &angular_step);
      /* write some data to a disk file. A header of 100 items will be written 
      first. I doubt if I will ever need 100 items, but its better to have too 
      many, rather than too few. */
      if(errno)
         fprintf(stderr, "Errno = %d in yagi.c 2\n", errno);
      ofp=fopen(output_filename,"wb");
      errno=0;
      write_header_to_disk(ofp, elements, driven, parasitic,min_frequency, max_frequency, frequency, step_frequency, angular_step);  
      /* write the x and y coordinates of the centre of the element, and the
      length in m */

      write_coordinates_of_elements_to_disk(ofp, driven, parasitic, driven_data, parasitic_data); 
      /* Now we have all the data in memory about the elements, we can calculate
      the impedance matrices. This will be a square matrix, of size NxN, where
      N is the number of elements in the yagi. The components on the diaganol
      Z11, Z22, Z33 etc, will be the self impedance of each element. The
      off-diagonal components, Zij (Z12, Z13 etc are the mutual impedance between
      elements. */
      current=FCOMPLEXvector(1L,(long) elements); 
      voltage=FCOMPLEXvector(1L,(long) elements);


      v=dvector(1L,2L*elements); 

      for(frequency=min_frequency; frequency <=max_frequency; frequency+=step_frequency)
      {
            fill_z_matrix(frequency,driven,parasitic,driven_data,parasitic_data, z);
            /* Now we must fill the V vector, which gives the voltage at the centre of
            each driven element. Since we know the magnitude of the voltage and phase,
            we can calculate the real and imaginary components of voltage. NB */
            fill_v_vector(driven, parasitic, driven_data, v);
            /* We now have the voltage vector V, and the impedance matrix Z. All
            we need to do now is solve a set of NxN equations, where N is the
            number of elements, to get N values of current, which we can put in
            the I matrix.
            Unfortunately, this is not trivual, by any standards!!!  It is
            complicated by the fact the the N equations are all complex. 
            The method used us a brute-force approach mentioned by Press
            in their 2nd Edition of the Numerical Recipes in C book. Here we
            put the data in the form:
      
            |Zr -Zc|   |Ic|   = |Vr|   which we will call A.x=b
            |Zc  Zr|   |Ic|     |Vc| 
            
            so the Z data goes now in a 2Nx2N matrix. This is a bit wasteful of space
            and time, but it will do here. */
            /* Read impedance data from 'z' to matrix A */
            for(i=1;i<=elements;++i)
            {
                  for(j=1;j<=elements;++j)
                  {
                        A[i][j]=z[i][2*j-1];   /* real data, top left corner of Z */
                        A[i+elements][j+elements]=z[i][2*j-1]; /* Bot right, real data */
                        A[i][j+elements]=-z[i][2*j]; /* imaginary, top right */
                        A[i+elements][j]=z[i][2*j]; /* bottom left, imaginary */
                  }
            }
            /* The following function prints to stantard output the z matrix of the
            antenna. It is really only used during debugging, so it can be commented
            out normally. */
            if(pflg)
                  print_z_matrix(frequency,elements,z);    
            /* read voltage data from v to b */
            for(i=1;i<=elements;++i)
            {
                  b[i]= v[2*i-1];                                                 /* real data */
                  b[i+elements]=v[2*i];                           /* complex data */
                  voltage[i].r=v[2*i-1];
                  voltage[i].i=v[2*i];
            /*
                  if(isnand(voltage[i].r) || isnand (voltage[i].i) || isnand(b[i]) || isnand(b[i+elements]) ||
                  (voltage[i].r) < 0.001)
                  {
                     fprintf(stderr, "b[%d]=%f b[%d]=%f %f %f\n", i, b[i], i+elements, b[i+elements], voltage[i].r, voltage[i].i);
                  }
            */
            }
            /* we now write the voltage data to disk, but only once, *not* at
            every frequency, since its fixed */
            if(frequency==min_frequency)
            {
                  fwrite((char *) &voltage[1], sizeof(struct FCOMPLEX), elements, ofp);
            }
            /* Perform a LU decompositon of A */
#ifdef DEBUG2
            for(test=1;test <= 2*elements; ++test)
            {
              for( =1;  <= 2*elements; ++)
                printf("A[%d][%d]=%f\n", test, , A[test][]);
            }
#endif
            ludcmp(A, elements*2, indx, &d); 
            /* We now have the voltages in b. After lubksb is run, we get the 
            currents in b */

            lubksb(A, 2*elements, indx, b); /* current's in b after lubksb has run*/
            /* Put currents an FCOMPLEX matrix current[element] 
            currents are stored in b as r,r,r,r ..... i,i,i,i */
            for(element_number=1;element_number<=elements;element_number++) 
            {
                  current[element_number].r=b[element_number];
                  current[element_number].i=b[element_number+elements];
            } 
            if(dflg)
                  display_antenna_currents(current,elements);
            fwrite((char *) &current[1], sizeof(struct FCOMPLEX), elements, ofp);
            real=b[1]/((b[1]*b[1])+(b[1+elements]*b[1+elements]));
            imaginary=-b[1+elements]/((b[1]*b[1])+(b[1+elements]*b[1+elements]));
         /* printf("f=%.2f MHz z = 1/(%f + i %f) = %f + i%lf \n", frequency/1e6, b[1], b[1+elements], real, imaginary);   */

            if(min_frequency < max_frequency && !sflg)
            {
                  printf("%s completed %5.1f%% f=%f MHz\n",argv[0],100*(frequency-min_frequency)/(max_frequency-min_frequency), frequency);
            }
      }
      if(errno)
         fprintf(stderr, "Errno = %d in yagi.c 3\n", errno);
      printf("\n");
      fclose(ofp);

      /* free strings */
      free_string(input_filename,0L,90L);
      free_string(line,0L,(long) MAX_LINE);
      free_string(output_filename,0L,90L);

      /* free vectors */
      free_dvector(x,1L, 2L*elements);
      free_dvector(b,1L, 2L*elements);
      free_dvector(v,1L, 2L*elements);
      free_ivector(indx,1,2L*elements);
      free_FCOMPLEXvector(current,1L,(long) elements); 
      free_FCOMPLEXvector(voltage,1L,(long) elements);
      /* free matrices */
      free_dmatrix(driven_data,1L,(long)(driven),1L,6L);
      // free_dmatrix(parasitic_data,1L,(long) parasitic,1L,4L); 
      free_dmatrix(z,1L,(long) elements, 1L, elements*2L);
      free_dmatrix(A,1L, 2L*elements, 1L, 2L*elements); 
      if(errno)
      {
         fprintf(stderr, "Errno = %d in yagi.c 4\n", errno);
         exit (errno);
      }
      else
            exit(0);
}

Generated by  Doxygen 1.6.0   Back to index