//NewGuide_State.java

import java.util.*;

public class NewGuide_State {
    //----------------------new variables for antenna
    public double theta_angle, theta_maximum, theta_minimum, 
		  phi_angle, phi_maximum, phi_minimum, 
		  point_distance, distance_minimum, distance_maximum, 
		  wavelength, frequency, wire_radius;
    public double Directivity, BeamWidth, TotalPower, TimeAveragePower;
    public double DipoleLength_lambda, DipoleLength_meters;
    public Complex Electric_Field, Magnetic_Field;
    public Complex ReflectionTE, ReflectionTM, ReflectionTE2, ReflectionTM2;
    public Double ReflectionMag;
    public double Tower1_height, Tower2_height;
    public double Tower_max, Tower_min;
    public double DirectPath, DirectPath_wavelength, BounceBP, BouncePD;
    public double BounceBP_wavelength, BouncePD_wavelength;
    public double BounceAP, BouncePC;
    public double BounceAP_wavelength, BouncePC_wavelength;
    public double PathDifference, PathDifference_meters, PathDifference_wavelength, 
            PhaseDifference, Approximation, fase;
    public double PhaseDifferencePath, PhaseDifferenceReflection;
    public double bounce_angle, bounce_angle_rad, max, Max_Grazing_Angle;
    public double dist, dist_min, dist_max;
    public boolean IsTE, onlyDirect;
    public double sigma, epsilon_r;
    public int GroundChoose;
    //--------------------------------------------------------------------------
    public double a, a_meters; // dimensions
    public double a_minimum, a_maximum, DeltaLength;
    public double mu_r, mu_r0, lambda_0, lambda; 
    public double taper_parameter, Maximus;
    
    public double DeltaZ;
    public double frequency_minimum, frequency_maximum, Deltaf;
    public double angular_frequency, beta_total;
    public double beta;
    public double phase_velocity, phase_velocity0, light_velocity, light_velocity0;
    
    public Complex ZeroC;
    public double distance;
    public double medium_impedance;
    public double DVector[], RVector[], RinVector[], XVector[], XinVector[], leng[], DVMax, RVMax, RinVMax, XVMax, XVMin, XinVMax, XinVMin;
    
    public double x[], x_polar[];
    public Complex etheta[];
    public double ethetaM[], ethetaM2[], ethetaM_dB[];
    public double EfieldMax, ethetaS, ethetaS2, ethetaS_dB; 
    public double area_front;
    public double zpos[], zposLarge[], zLambda[], theta, phi, zposMax;
    
    public int Numpoints, Numx, Nsections, NsectionsLarge;
    public double dz, angle_factor;
    public boolean inputpanelflag, Is_theta, Is_radius, IsPolar, Is_dB;
    public boolean microFlag2, microFlag3;
    
    public int Nmaxint, ShowWhich;
    
    public boolean IsDirectivity, IsRad, IsRin, IsXrad, IsXin;
    
    public static final double epsilon0 = 8.841941286E-12; // approximate value
    public static final double mu0 = 4.0*Math.PI*1.0e-7;//1.25663706144E-6; //Units H/m
    public Complex jay = new Complex(0.0,1.0);
    public Complex minusjay = new Complex(0.0,-1.0);
    public Complex One = new Complex(1.0,0.0);
    public Complex minusOne = new Complex(-1.0,0.0);
    
    public boolean LicenseExpired;
    
    public int this_month, today_week, this_year, this_hour, this_minute, today_month, 
	       today_year,this_zone, saving_time;
	       
    GregorianCalendar Greg = new GregorianCalendar();

    
    public NewGuide_State(){
    
	this_month = Greg.get(Calendar.MONTH);
	today_week = Greg.get(Calendar.DAY_OF_WEEK);
	this_year = Greg.get(Calendar.YEAR);
	this_hour = Greg.get(Calendar.HOUR_OF_DAY);
	this_minute = Greg.get(Calendar.MINUTE);
	today_month = Greg.get(Calendar.DAY_OF_MONTH);
	today_year = Greg.get(Calendar.DAY_OF_YEAR);
	this_zone = Greg.get(Calendar.ZONE_OFFSET);
	saving_time = Greg.get(Calendar.DST_OFFSET);
	
	wire_radius = 1.0e-5; //radius of dipole wire (wavelengths) 
	
	// Tables for Sin and Cos Integrals
	Nmaxint = 1001;
        ShowWhich = 1;
        max = 0;
        ReflectionTE = new Complex(-1.0,0.0);
        ReflectionTM = new Complex(1.0,0.0);
        ReflectionMag = 1.0;
	//new initialization for antenna
	    
            taper_parameter = 0.5;
            
	    theta_angle = 135.0;//degrees
	    theta_maximum = 180.0;
	    theta_minimum = 0.0;
	    
	    phi_maximum = 360.0;
	    phi_minimum = 0.0;
	    point_distance = 10.0; //wavelength
	    distance_minimum = 0;//      "
	    distance_maximum = 1000;//   "
	    
            //======================================
            medium_impedance = Math.sqrt(mu0*mu_r/(epsilon0*epsilon_r));
            frequency = 1.0E9;
            frequency_minimum = 30.0E6;
            frequency_maximum = 3.0E9;
            angular_frequency = 2.0 * Math.PI * frequency;

            lambda_0 = 1.0/Math.sqrt(epsilon0*mu0)/frequency;
            lambda = lambda_0;
            ZeroC = new Complex(0.0,0.0);
            light_velocity = 1.0/(Math.sqrt(epsilon0*epsilon_r*mu0*mu_r));
            light_velocity0 = 1.0/(Math.sqrt(epsilon0*mu0));
            wavelength = 1.0/Math.sqrt(epsilon0 * mu0)/frequency;//meters
            
            dist = 5000.0; //meters
	    dist_min = 500.0;//      "
	    dist_max = 10000.0;//   "
            
            frequency = 1.0E9;//Hertz
	    Tower1_height = 100.0; // meters
            Tower2_height = 50.0; // meters
        
            Tower_min = 0.0;// meters
            Tower_max = 500.0;// meters
            
            DirectPath = Math.sqrt(dist*dist + Math.pow(Tower2_height-Tower1_height,2));
            DirectPath_wavelength = DirectPath/wavelength;
            BounceBP = Tower1_height*dist/(Tower1_height+Tower2_height);
            BouncePD = dist - BounceBP;
            
            BounceBP_wavelength = BounceBP/wavelength;
            BouncePD_wavelength = BouncePD/wavelength;
            
            BounceAP = Math.sqrt(BounceBP*BounceBP + Math.pow(Tower1_height,2));
            BouncePC = Math.sqrt(BouncePD*BouncePD + Math.pow(Tower2_height,2));
            BounceAP_wavelength = BounceAP/wavelength;
            BouncePC_wavelength = BouncePC/wavelength;
            
            PathDifference_meters = BounceAP+BouncePC-DirectPath;
            PathDifference_wavelength = PathDifference/wavelength;
            bounce_angle_rad = Math.atan((double)Tower1_height/(double)BounceBP);
            bounce_angle = bounce_angle_rad*180/Math.PI;
            
            double Try;
            Try = PathDifference_wavelength;
            while(Try > 1.0){Try -= 1.0;}
        
            PhaseDifferencePath = Try * 360; // degrees
            //System.out.println(" Phase-Path = "+PhaseDifferencePath);
        
            sigma = 1.0E130;
            GroundChoose = 4;
        
            //----------------------------------------------------------------------
            double Find;
            if(IsTE){
                Find = PathDifference_wavelength + 0.5;
            }
            else{
                Find = PathDifference_wavelength;
            }
        
            while(Find > 1.0){
                Find -= 1.0; 
            }
            //----------------------------------------------------------------------    
            
            PhaseDifference = - Find*360.0; // degrees
            phi_angle = - Find*2.0*Math.PI; // radians
            Max_Grazing_Angle = 0.0;
            
            Approximation = 2*Tower1_height*Tower2_height/dist;
            IsTE = true; onlyDirect = false;
            
            //ReflectionTE2 =  EMF.getOblique_Reflection_Coefficient((Math.PI*0.5 - bounce_angle_rad), 1.0, epsilon_r,
		//1.0, 1.0, 0.0, sigma, angular_frequency, false);
            
            //ReflectionTM2 = EMF.getOblique_Reflection_CoefficientH((Math.PI*0.5 - bounce_angle_rad), 1.0, epsilon_r,
		//1.0, 1.0, 0.0, sigma, angular_frequency, true);
            
            //System.out.print(ReflectionTE2);
            
            //=======================================
            area_front = distance * distance * 4.0 * Math.PI; // Area of spherical front
            
	    Numpoints = 1001;
            Numx = 1001;
            
	    a = 2.0; // wavelengths of antenna length
	    DipoleLength_lambda = a;//lambda
            
	    a_minimum = 2.0;
	    a_maximum = 12.0;// wavelengths
	    DeltaLength = (a_maximum - a_minimum)/(Numpoints - 1);
	    
	    DipoleLength_meters = DipoleLength_lambda * wavelength;
	    Directivity = 0.0;
	    TotalPower = 0.0;
	    
	    DVector = new double[Numx];
	    leng = new double[Numx];
	    
	//------------------------------------------------------------------
	
	x_polar = new double[361];
	for(int i=0; i<361; i++){
	    x_polar[i] = i;
	}

	epsilon_r = 1.0;
	mu_r = 1.0;
	
	angle_factor = Math.PI/180.0;
	
	medium_impedance = Math.sqrt(mu0*mu_r/(epsilon0*epsilon_r));
	frequency = 1.0E9;
	frequency_minimum = 30.0E6;
	frequency_maximum = 3.0E9;
	angular_frequency = 2.0 * Math.PI * frequency;
	
	lambda_0 = 1.0/Math.sqrt(epsilon0*mu0)/frequency;
	lambda = lambda_0;
	ZeroC = new Complex(0.0,0.0);
	light_velocity = 1.0/(Math.sqrt(epsilon0*epsilon_r*mu0*mu_r));
	light_velocity0 = 1.0/(Math.sqrt(epsilon0*mu0));
	
	Nsections = 40;
	NsectionsLarge = 501;
	
	phi = Math.PI/2.0;
	
	inputpanelflag = true;
	Deltaf = (frequency_maximum - frequency_minimum)/(Numpoints - 1);
	
	Is_theta = true;
	Is_radius = false;
	IsPolar = false;
	IsDirectivity = true;
	IsRad = false;
	IsRin = false;
	IsXrad = false;
	IsXin = false;
	Is_dB = true;
	microFlag2 = false; microFlag3 = false; 
        double pr = 0;
        
	//scan_coefficients();
	//ignition();     
        get_maxangle();
        
 }
    
    public void ignition(){
	int Lratio;
	angular_frequency = 2.0 * Math.PI * frequency;
	
	light_velocity0 = 1.0/(Math.sqrt(epsilon0*mu0));
        light_velocity = 1.0/(Math.sqrt(epsilon0*mu0));
        
        lambda_0 = light_velocity0/frequency;
	lambda = light_velocity/frequency;
	
	beta = 2.0*Math.PI/lambda;
	medium_impedance = Math.sqrt(mu0*mu_r/(epsilon0*epsilon_r));
	
        if(GroundChoose == 1){
            epsilon_r = 1.0;
            sigma = 1.0E130;
        }
        else if(GroundChoose == 2){
            epsilon_r = 5.0;
            sigma = 1.0E-3;
        }
        else if(GroundChoose == 3){
            epsilon_r = 15.0;
            sigma = 5.0E-3;
        }
        else if(GroundChoose == 4){
            epsilon_r = 30.0;
            sigma = 2.0E-2;
        }
        else if(GroundChoose == 5){
            epsilon_r = 81.0;
            sigma = 1.0e-2;
        }
        else if(GroundChoose == 6){
            epsilon_r = 81.0;
            sigma = 5.0;
        }
        
        wavelength = 1.0/Math.sqrt(epsilon0 * mu0)/frequency; //meters
        DirectPath = Math.sqrt(dist*dist + Math.pow(Tower2_height-Tower1_height,2));
        DirectPath_wavelength = DirectPath/wavelength;
        BounceBP = Tower1_height*dist/(Tower1_height+Tower2_height);
        BouncePD = dist - BounceBP;

        BounceBP_wavelength = BounceBP/wavelength;
        BouncePD_wavelength = BouncePD/wavelength;

        BounceAP = Math.sqrt(BounceBP*BounceBP + Math.pow(Tower1_height,2));
        BouncePC = Math.sqrt(BouncePD*BouncePD + Math.pow(Tower2_height,2));
        BounceAP_wavelength = BounceAP/wavelength;
        BouncePC_wavelength = BouncePC/wavelength;

        PathDifference_meters = BounceAP+BouncePC-DirectPath;
        PathDifference_wavelength = PathDifference_meters/wavelength;
        
        if(BounceBP<BouncePD){
            bounce_angle = Math.atan((double)Tower2_height/(double)BouncePD)*180.0/Math.PI;
        }
        else{
            bounce_angle = Math.atan((double)Tower1_height/(double)BounceBP)*180.0/Math.PI;
        }
        bounce_angle_rad = bounce_angle * Math.PI/180.0;
        
        double Try;
        Try = PathDifference_wavelength;
        while(Try > 1.0){
            Try -= 1.0; 
        }
        PhaseDifferencePath = - Try * 360; // degrees
        //----------------------------------------------------------------------
        // Find the Reflection coefficients
        double coseno, seno, sigmafactor;
        sigmafactor = 0.0;
        
        if(GroundChoose == 1){
            ReflectionTE = new Complex(-1.0,0.0);
            ReflectionTM = new Complex(1.0,0.0);
        }
        else{
            
            Complex rootfactor, numeratorTE, denominatorTE, numeratorTM, denominatorTM;
            coseno = Math.cos(bounce_angle_rad);
            seno = Math.sin(bounce_angle_rad);
            sigmafactor = - sigma * 0.5/(Math.PI*frequency*epsilon0);
            
            rootfactor = Complex.Pow(new Complex(epsilon_r-coseno*coseno,sigmafactor), 0.5);
            numeratorTE = Complex.Subtract(new Complex(seno,0.0), rootfactor);
            denominatorTE = Complex.Add(new Complex(seno,0.0), rootfactor);
            
            numeratorTM = Complex.Subtract(new Complex(seno*epsilon_r,seno*sigmafactor), rootfactor);
            denominatorTM = Complex.Add(new Complex(seno*epsilon_r,seno*sigmafactor), rootfactor);
                     
            ReflectionTE = Complex.Divide(numeratorTE, denominatorTE);
            ReflectionTM = Complex.Divide(numeratorTM, denominatorTM);
            
           // ReflectionTE2 =  EMF.getOblique_Reflection_Coefficient((Math.PI*0.5 - bounce_angle_rad), 1.0, epsilon_r,
	   //	1.0, 1.0, 0.0, sigma, angular_frequency, false);
            
           // ReflectionTM2 = EMF.getOblique_Reflection_CoefficientH((Math.PI*0.5 - bounce_angle_rad), 1.0, epsilon_r,
	   //	1.0, 1.0, 0.0, sigma, angular_frequency, true);
            
            //System.out.print("  GammaTE2 = "+ReflectionTE2+"   GammaTM2 = "+ReflectionTM);

        }
        
        //----------------------------------------------------------------------
        double Find;
        Find = 0.0;
        fase = 0.0;
        
        if(GroundChoose == 1){
            if(IsTE){
                Find = PathDifference_wavelength + 0.5;
                ReflectionMag = 1.0;
            }
            else{
                Find = PathDifference_wavelength;
                ReflectionMag = 1.0;
            }
        }
        else{
            if(IsTE){
                fase = Complex.Arg2(ReflectionTE);
            
                Find = PathDifference_wavelength - fase/(2.0*Math.PI);
                ReflectionMag = Complex.Magnitude(ReflectionTE);
                //System.out.println("TE  - Dphi = "+ PathDifference_wavelength+"   fase = "+360*fase/(2.0*Math.PI));
            }
            else{
                fase = Complex.Arg2(ReflectionTM);
                Find = PathDifference_wavelength - fase/(2.0*Math.PI);
                ReflectionMag = Complex.Magnitude(ReflectionTM);
                //System.out.println("TM  - Dphi = "+ PathDifference_wavelength+"   fase = "+360*fase/(2.0*Math.PI));
            }
        }
        
        while(Find > 1.0){
            Find -= 1.0; 
        }
        //----------------------------------------------------------------------    
            
        PhaseDifference = - Find*360.0; // degrees
        phi_angle = - Find*2.0*Math.PI; // radians
        //System.out.println("phase = "+fase+"  phi = "+phi_angle);
        
        Approximation = 2*Tower1_height*Tower2_height/dist;
        
        if(Tower1_height == 0.0 || Tower2_height == 0.0){
            onlyDirect = true;
        }
        else{
            onlyDirect = false;
        }
         
        get_maxangle();
        
        Max_Grazing_Angle = MaestroA.rounder(90.0 - max/100.0,4);
        
        //System.out.print("max angle = "+Max_Grazing_Angle);
        //System.out.println("GammaTE = "+ReflectionTE+"  GammaTM = "+ReflectionTM);
        //System.out.println("sigma/omega*epsilon"+sigmafactor+"    ReflectionMag = "+ReflectionMag);
        
        //----------------------------------------------------------------------
        DipoleLength_lambda = a;//lambda
	DipoleLength_meters = DipoleLength_lambda * wavelength;
	Lratio = (int)(100*DipoleLength_meters/wavelength)+1;
	if(Lratio < 51){Lratio = 51;}
	Nsections = Lratio;
	//NsectionsLarge = Lratio*10+1;
	NsectionsLarge = 501;
	DeltaZ = DipoleLength_meters / Nsections;
	DeltaLength = (a_maximum - a_minimum)/(Numpoints - 1);
	
	distance = point_distance*wavelength;//meters
	area_front = distance * distance * 4.0 * Math.PI;
	
        etheta = new Complex[361];
	ethetaM = new double[361];
	ethetaM2 = new double[361];
        ethetaM_dB = new double[361];
        
        zpos = new double[Nsections];
	dz = DipoleLength_meters / Nsections;
	
	angular_frequency = 2.0 * Math.PI * frequency;
	
	light_velocity0 = 1.0/(Math.sqrt(epsilon0*mu0));
        light_velocity = 1.0/(Math.sqrt(epsilon0*mu0));
	lambda_0 = light_velocity0/frequency;
	lambda = light_velocity/frequency;
	
	beta = 2.0*Math.PI/lambda;
	medium_impedance = Math.sqrt(mu0*mu_r/(epsilon0*epsilon_r));
		
	// Obtain Directivity, Radiation Resistance and Time-Average Power  (Far Field assumption)
	double Func[], Func2[], Func2Max, angle, betafactor, FuncInt, Delta,theta_angle_radians;
	int imax, Nmax;
	Nmax = 1800;
	Func   = new double[Nmax+1];
	Func2   = new double[Nmax+1];
	Func2Max = 0.0;
	imax = 0;
	Delta = 0.5*Math.PI/Nmax;
	betafactor = beta*DipoleLength_meters*0.5;
	theta_angle_radians = theta_angle*Math.PI/180.0;
	Func[0] = 0.0;
	Func2[0] = 0.0;
        
        //scanField();
        //double Maxim = MaestroA.getMax(ethetaM);
        // scan theta ----------------------------------------------------------
        double sum_thetaS;
        double factorAS;
        double xmaxS, DeltaxS, Maxx;
        xmaxS = DipoleLength_meters * 0.5;
	DeltaxS = xmaxS / Nsections;
	Maxx = 0.0;
        factorAS = 0.0;
        for(int j=0;j<Nsections+1;j++){

            factorAS = (1.0-DeltaxS*j/xmaxS*taper_parameter)*DeltaxS*j
                      *MaestroA.J0(DeltaxS*j*beta * Math.sin(0.0));    
            Maxx += factorAS;
        }
        Maxx = Maxx/xmaxS;
        factorAS = 0.0;
        sum_thetaS = 0.0;
        for(int j=0;j<Nsections+1;j++){    
            factorAS = (1.0-DeltaxS*j/xmaxS*taper_parameter)*DeltaxS*j
                      *MaestroA.J0(DeltaxS*j*beta * Math.sin(theta_angle_radians-(Math.PI*0.5)));    
            sum_thetaS += factorAS;
        }
        sum_thetaS = sum_thetaS/xmaxS;
        ethetaS = Math.abs(sum_thetaS)/Maxx;
        ethetaS2 = Math.pow(ethetaS,2.0);
        ethetaS_dB = 10.0*Math.log10(ethetaS2);
        //----------------------------------------------------------------------
	//scanField();
        //scan_coefficients();
}

public void scanField(){
    double xstart, xmax, Deltax, sum_theta, sum_dir;
    double theta_deg, theta_rad, factorA;
    int bw;
	xstart = 0.0;
        bw = 0;
        xmax = DipoleLength_meters * 0.5;
        Deltax = xmax / Nsections;
	theta = theta_angle*angle_factor; 
	factorA = 0.0;
	Maximus = 0.0;
        for(int i=0;i<91;i++){
            theta_deg = (double)i;
            theta_rad = theta_deg * angle_factor;
            // simple integration 
            sum_theta = 0.0;
            double factor = 0.0;
            for(int j=0;j<Nsections+1;j++){
                factor = Deltax*j*beta * Math.sin(theta_rad);
                factorA = (1.0-Deltax*j/xmax*taper_parameter)*Deltax*j
                         *MaestroA.J0(factor); 
                sum_theta += factorA;
            }
            
            sum_theta = (Math.abs(sum_theta/xmax));
            ethetaM[i+90] = Math.abs(sum_theta);
            ethetaM[90-i] = Math.abs(sum_theta);
        }
        Maximus = MaestroA.getMax(ethetaM);
        
        double temp = 0.0;
        bw = 0;
        for(int i=0;i<181;i++){
            temp = ethetaM[i]/Maximus;
            ethetaM[i] = (temp);
            ethetaM2[i] = temp*temp;
            ethetaM_dB[i] = 10.0*Math.log10(temp*temp);
            if(ethetaM_dB[i] < -100.0){ethetaM_dB[i] = -100.0;}
        }
        
        // Calculate Directivity
        sum_dir = 0.0;
        for(int j=0;j<181+1;j++){
            //sum_dir += ethetaM2[j]*Math.sin(j*Math.PI/180.0f)*Math.PI/180.0f;// * Math.sin((double)j*angle_factor));  
            sum_dir += ethetaM2[j];//*Math.PI/180.0f;// * Math.sin((double)j*angle_factor));  
        }
        Directivity = 4*Math.PI/(Math.abs(sum_dir*angle_factor)*Math.abs(sum_dir*angle_factor));
        //Directivity = 2.0/Math.abs(sum_theta*angle_factor);
        
        // Calculate numerically beamwidth
        bw = 0;
        for(int i=90;i<181;i++){
            if(ethetaM_dB[i] < -3.0){bw = i-90; break;}
        }
        
        // interpolate for beamwidth
        double DeltaE, DeltaBW, extra;
        DeltaE = Math.abs(ethetaM2[bw+90-1] - ethetaM2[bw+90]);
        DeltaBW = Math.abs(0.5 - ethetaM2[bw+90-1]);
        extra = DeltaBW/DeltaE;
        BeamWidth = 2.0*((double)(bw-1) + extra);
        
}        
            
        

//---------------------------


    public void scan3Dfield(){
        
    }
    public void scan_reactance(){
	
    }


    //Gets array for plot of Field as function of frequency    
    public void scan_coefficients(){
	
	DVMax = 0.0;
	for(int iL=0;iL<Numx;iL++){
            leng[iL] = iL;
            DVector[iL] = 1.0 - taper_parameter*iL/(Numx-1);	
	}
    }
    
    public void get_maxangle(){
	    max = EMF.get_Trans_max(1.0,epsilon_r,1.0,1.0,0.0,sigma,angular_frequency,true);
    }
}