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

gain.c

#ifdef HAVE_STDLIB_H
#include <stdlib.h>
#endif
#include <stdio.h>
#include <math.h>
#include <malloc.h>
#include <errno.h>
#include "yagi.h"
#define TINY 1e-10


/* This function finds the gain both in the E plane (xz plane) and the
H plane (xy plane) at angle (theta, phi). The method used is as described
on page 1-12 of 'Yagi Antenna Design' by Dr. Lawson , ARRL */

void gain(double theta, double phi, double pin, double F, struct element_data *coordinates, struct FCOMPLEX *current, int elements, double *gain_E_plane, double *gain_H_plane, double actual_frequency, double design_frequency) 
{
      int i;
      double *r_E, *r_H, *g_E, *g_H,integer_bit;
      double length, x, y, lamda_design, lamda, tmp;
      struct FCOMPLEX temp_E, temp_H, e_gain, h_gain;
      /* need to allocate space for FCOMPLEX types. Since there is no Numerical 
      Recipes routine, I'll use the standard method. Since I've always used
      elements positions 1 to N, I'll just wast the extra location */
      r_E=dvector(1L,(long) elements); 
      r_H=dvector(1L,(long) elements);
      g_E=dvector(1L,(long) elements);
      g_H=dvector(1L,(long) elements);
      e_gain.r=0;   /* set to zero. Necessary as we sum into these */
      e_gain.i=0;
      h_gain.r=0;
      h_gain.i=0;
      /* convert theta and thi to radians from degrees */
      theta=theta*M_PI/180;
      phi=phi*M_PI/180;
      lamda_design=3e8/design_frequency;
      lamda=3e8/actual_frequency;
      for(i=1;i<=elements;++i) 
      {
            length=coordinates[i].length/lamda;
            x=coordinates[i].x/lamda_design;
            y=coordinates[i].y/lamda_design;

            /* for E -plane */
            r_E[i]=sin(theta)*x;
            if(fabs(theta) < TINY)    /* avoid division by zero if theta=0 */
                  g_E[i]=0;
            else
                  g_E[i]=(cos(M_PI*length*cos(theta))-cos(M_PI*length))/sin(theta);

            /* for H -plane */
            r_H[i]=x*cos(phi)+y*sin(phi);
            g_H[i]=1-cos(M_PI*length);
            /* printf("g_H[%d]=%.16lf \n", i, g_H[i]);             */
            temp_E.r=0;
            temp_E.i=2*M_PI*r_E[i]*F;
            temp_E=E_to_complex_power(temp_E); /* exp(j 2 pi r[i] F) */

            temp_H.r=0;
            temp_H.i=2*M_PI*r_H[i]*F;
            temp_H=E_to_complex_power(temp_H);
      /* printf("element %d temp_H.r=%f temp_H.i=%f\n",i, temp_H.r, temp_H.i); */
            /* get element currents */
            temp_E=Cmul(temp_E,current[i]);
            temp_H=Cmul(temp_H,current[i]);
            /* printf("element %d: current = %.10lf i%.10lf = %.10lf (mag) at %f degrees\n", i,current[i].r, current[i].i, sqrt(current[i].r*current[i].r+current[i].i*current[i].i), 180*atan2(current[i].i,current[i].r)/M_PI);         */
      /* printf("element %d I temp_H.r=%.10lf temp_H.i=%.10lf\n",i, temp_H.r, temp_H.i); */
            e_gain=Cadd(e_gain,RCmul(g_E[i],temp_E));
            h_gain=Cadd(h_gain,RCmul(g_H[i],temp_H));
            /* printf("h_gain=%.16lf + i %.16lf \n", h_gain.r, h_gain.i); */
      }
      /* there will be a divide  by zero in calculating
      equ 1.28, of Lawsons book, (see page 1-12)
      at theta=0,180,360 degree, infact anywhere where 
      sin(theta) is zero */

      if( modf(theta/M_PI,&integer_bit) < TINY || fabs(theta-M_PI)<TINY)
      {
            *gain_E_plane=-150.0;   /* a very very small number */
      }
      else
      {
            tmp=(Cabs(e_gain)*Cabs(e_gain)*60/pin);  
            *gain_E_plane=10*log10(tmp);  
#ifdef DEBUG
      if(errno)
      {
            fprintf(stderr,"Errno =%d in gain() of gain.c after E\n", errno);
            printf("tmp=%f pin=%f theta=%f gainE=%f\n", tmp, pin,theta,*gain_E_plane);
            exit(1);
      }
#endif
      }
      tmp=(Cabs(h_gain)*Cabs(h_gain)*60/pin);  
      *gain_H_plane=10*log10(tmp);  

#ifdef DEBUG
      if(errno)
      {
            fprintf(stderr,"Errno =%d in gain() of gain.c after H\n", errno);
            printf("tmp=%f pin=%f \n", tmp, pin);
            exit(1);
      }
#endif
      free_dvector(r_E,1L,(long) elements); 
      free_dvector(r_H,1L,(long) elements); 
      free_dvector(g_E,1L,(long) elements); 
      free_dvector(g_H,1L,(long) elements); 

#ifdef DEBUG
      if(errno)
      {
            fprintf(stderr,"Errno =%d in gain() of gain.c\n", errno);
            exit(1);
      }
#endif
}

struct FCOMPLEX E_to_complex_power(struct FCOMPLEX x)
{
      struct FCOMPLEX a, answer;
      double real;

      /* First real part */
      real=exp(x.r);
      /* now the imaginary part */
      /* Exp(i theta) = cos(theta) + i sin(theta) */
      a.r=cos(x.i);
      a.i=sin(x.i);
      answer=RCmul(real,a);
      return(answer);
}
#undef TINY

Generated by  Doxygen 1.6.0   Back to index