process.c

FDTD 模拟 天线 包含 近场－远场 外推程序（时谐场）

/* ------------------------------------------ 

     The Result Processing ( DFT )

*/

# include "process.h"


/* Calculating the Voltage between two Points (X_Location,Y_Location,Z_start)
  
   and (X_Location,Y_Location,Z_end) 

   by Discrete Fourier Transform 

*/

void DFT_Volt(double **** Ez,int time, int nX_Location, int nY_Location,
              int nZ_start, int nZ_end, double *V_R, double *V_I)
{

  int i,k;

  double V_tmp_1;

  V_tmp_1=0.0;  

  for(k=nZ_start; k<nZ_end; k++)
   {
         
    V_tmp_1+=Ez[1][nX_Location][nY_Location][k]*dhz;

    }        

   for(i=0;i<N_NT; i++)
    {
       
      V_R[i]+=V_tmp_1*cos(2.0*Pi*i*time/N_NT);
      V_I[i]+=-V_tmp_1*sin(2.0*Pi*i*time/N_NT);


    }  
       

}     


/*

    Exchange the Value of Field at the last time E[0]
and current time E[1]

*/
  

void Exchange_Value(double ****Ex, double ****Ey, double ****Ez, double ****Hx,
                    double ****Hy, double ****Hz)
{

  int i,j,k;

  int NExx,NExy,NExz,NHxx,NHxy,NHxz;
  int NEyx,NEyy,NEyz,NHyx,NHyy,NHyz;
  int NEzx,NEzy,NEzz,NHzx,NHzy,NHzz;

  double temp;

  NExx=Nx;   NExy=Ny+1; NExz=Nz+1;
  NEyx=Nx+1; NEyy=Ny;   NEyz=Nz+1;
  NEzx=Nx+1; NEzy=Ny+1; NEzz=Nz;
  NHxx=Nx+1; NHxy=Ny;   NHxz=Nz;
  NHyx=Nx;   NHyy=Ny+1; NHyz=Nz;
  NHzx=Nx;   NHzy=Ny;   NHzz=Nz+1;

/* ------------- Ex --------------- */

  for(i=0;i<NExx; i++)
    {
      for(j=0;j<NExy; j++)
       {
         for(k=0;k<NExz; k++)
          {
            temp=Ex[0][i][j][k];
            Ex[0][i][j][k]=Ex[1][i][j][k];
            Ex[1][i][j][k]=temp;
          }

       }
   }

/* --------------- Ey ----------------- */

  for(i=0;i<NEyx; i++)
    {
      for(j=0;j<NEyy; j++)
       {
         for(k=0;k<NEyz; k++)
          {
            temp=Ey[0][i][j][k];
            Ey[0][i][j][k]=Ey[1][i][j][k];
            Ey[1][i][j][k]=temp;
          }

       }
   }

/* ----------------- Ez -------------- */
          
  for(i=0;i<NEzx; i++)
    {
      for(j=0;j<NEzy; j++)
       {
         for(k=0;k<NEzz; k++)
          {
            temp=Ez[0][i][j][k];
            Ez[0][i][j][k]=Ez[1][i][j][k];
            Ez[1][i][j][k]=temp;
          }

       }
   }

/* ------------------------ Hx -------------------- */

   for(i=0;i<NHxx; i++)
    {
      for(j=0;j<NHxy; j++)
       {
         for(k=0;k<NHxz; k++)
          {
            temp=Hx[0][i][j][k];
            Hx[0][i][j][k]=Hx[1][i][j][k];
            Hx[1][i][j][k]=temp;
          }

       }
   } 

/*  ------------------------- Hy -------------------- */

   for(i=0;i<NHyx; i++)
    {
      for(j=0;j<NHyy; j++)
       {
         for(k=0;k<NHyz; k++)
          {
            temp=Hy[0][i][j][k];
            Hy[0][i][j][k]=Hy[1][i][j][k];
            Hy[1][i][j][k]=temp;
          }

       }
   } 

/*   ------------------------ Hz ------------------------- */

   for(i=0;i<NHzx; i++)
    {
      for(j=0;j<NHzy; j++)
       {
         for(k=0;k<NHzz; k++)
          {
            temp=Hz[0][i][j][k];
            Hz[0][i][j][k]=Hz[1][i][j][k];
            Hz[1][i][j][k]=temp;
          }

       }
   } 


} 

void Cal_2D_Para(MODULE *pModule,char *pFileName, double *V_1_R, double *V_1_I,
                  double *V_2_R, double *V_2_I)
{
  int i;  
  double f,beta,Eref; 
  FILE *fp1,*fp2;

  fp1=fopen("Imp.dat","w+");
  fp2=fopen(pFileName,"w+");

/* --------- Record the Incidence Voltage ------- */

   for(i=0;i<N_NT; i++)
    {
     
      fprintf(fp2,"%lf  %lf \n",V_1_R[i],V_1_I[i]);
   
    }

    fclose(fp2);     

/* ===============================================  */


  fprintf(fp1,"====================================\n");
  fprintf(fp1," Frequency(GHz)      Equivalent Permittivity   \n");

  for(i=0;i<N_NT/2; i++)
   {
      f=i/(N_NT*dt);   
       
      beta=Phase_Cal(pModule,V_1_R[i],V_1_I[i],V_2_R[i],V_2_I[i]);

      Eref=beta*beta/(4.0*Pi*Pi*f*f*u0*e0);
     
      f=f/1.0e9;    /* Frequency (GHz) */  

      fprintf(fp1,"%lf %lf\n",f,Eref);


   }


   fclose(fp1);
   
}


  

   

double Phase_Cal(MODULE *pModule,double v1R,double v1i,double v2R,double v2i)
{

  double tmpr,tmpi,tmp;
  double dis;

  dis=N_2D_Reference*dhy;

  tmpr=(v1R*v2R+v1i*v2i)/(v2R*v2R+v2i*v2i);
  tmpi=(v1i*v2R-v1R*v2i)/(v2R*v2R+v2i*v2i);

  tmpr=(v2R*v1R+v2i*v1i)/(v1R*v1R+v1i*v1i);
  tmpi=(v2i*v1R-v2R*v1i)/(v1R*v1R+v1i*v1i);
  if(tmpi >0 )  
    {
       tmp=2.0*Pi+atan2(-tmpi,tmpr);
     }
   else
     {
        tmp=atan2(-tmpi,tmpr);
      }

   tmp=tmp/dis;  /* Phase Constant */

   return tmp;

}


void Cal_3D_Para(MODULE *pModule,char *pFileName, double *V_1_R, double *V_1_I)
{
  int i;  
  double f,S11; 
  double *V_in_R, *V_in_I;

  FILE *fp1,*fp2;

  fp1=fopen(pFileName,"r");

  V_in_R=dmalloc_1(N_NT); V_in_I=dmalloc_1(N_NT);

/* --------- Read the Incidence Voltage ------- */

   for(i=0;i<N_NT; i++)
    {
     
      fscanf(fp1,"%lf%lf \n",&V_in_R[i],&V_in_I[i]);
   
    }

    fclose(fp1);     

/* =============================================== 

     Calculate the S-Parameter and Record 

*/

  fp2=fopen("S.dat","w+");

  fprintf(fp2, "====================================\n");
  fprintf(fp2, " Frequency(GHz)   |S11| (dB)     \n");

  for(i=0;i<N_NT/2; i++)
   {
      f=i/(N_NT*dt);   
       
    S_Para(V_1_R[i],V_1_I[i],V_in_R[i],V_in_I[i],&S11);
     
   f=f/1.0e9;    /* Frequency (GHz) */  

   fprintf(fp2,"%lf %lf  %lf\n",f,S11);

   }

   fclose(fp2);

   dfree_1(V_in_R); dfree_1(V_in_I);
   
}  


void S_Para(double V_1_R,double V_1_I,double V_in_R,double V_in_I,double *pS11)
{

  double tmp,tmp_r,tmp_i;

  tmp=10.0*log10(V_in_R*V_in_R+V_in_I*V_in_I);

  tmp_r=V_1_R-V_in_R; tmp_i=V_1_I-V_in_I; 

  *pS11=10.0*log10(tmp_r*tmp_r+tmp_i*tmp_i)-tmp;


}