//-----------------------------------------------------------------------------
// 
// ORIGINAL AUTHOR: Steven Schultz, WB8WGY
// 
// DESCRIPTION:
//
//		Program which computes g(i) prototype values and their scaled network
//		values for broadband matching.
//
//		Part 1: Preliminary
//				Band Center Frequency and Slope Parameters are calculated.
//		Part 2: Fano's Broadband-Matching Limitations
//				Levy Matching
//				Newton - Raphson Solution
//				Load Decrement and number of elements (N) used to calculate 
//				parameters A and B which in turn are used to caculate 
//				SWR min and max and Return Loss min and max.
//		Part 3: Resistive Source with g(i) Prototype Values.
//				Load Decrement, A, B, and N parameters used to calculate
//				g(i) prototype values.
//		Part 4: Network Element Values Calculated.
//				g(i) prototype values are impedance or adittance scaled to
//				create final network values.
//
//		Basis of this software is Program B6-3.	
//		Newton - Raphson Solution then Levy Matching,
//		which is in Appendix B in book:
//
//			Circuit Design Using Personal Computers
//			Thomas R. Cuthbert, Jr.
//			1983 John Wiley and Sons 
//
//		BASIC program B6-3 is found on page 443 of this text.
//		The BASIC program is the basis for this C++ source file.  
//		A replica of this program is at the end of this file.
//
//-----------------------------------------------------------------------------
//
// Revision History:
//
// Rev: n January 30, 2009
//
//		Load Decrement and Load Q Limits revised.
//		Program comments and output text revised.
//		Added warnings for Load Decrement and Load Q Limits exceeded when
//		load parameters are used.
//		Provide output of date and time of program execution.
//
// Rev: m January 21, 2009
//
//		New input selections:
//		Three types of input made available to user:
//			1) Load Decrement obtained from Load Parameters
//			2) Load Decrement entered directly
//			3) Load Decrement obtained from Q and Fractional Bandwidth
//		The user then selects:
//			1) Series Load or Series Equivalent Load Parameters
//			2) Parallel Load or Parallel Equivalent Load Parameters
//
//		The R, L, and C elemental values are output for all types of input.
//
//		Added program file name, revision, and revision date to beginning of 
//		all program outputs: console, text, and csv.
//
// Rev: l January 3, 2009
//
//		Program calculates element values:
//			1) Load Reistance, G(0)
//			2) Elemental L and C for each G(1) through G(N)
//			3) Source Resistance, G(N+1)			
//
// Rev: k December 18, 2008
//
//		Program modified so one of four types of entry are allowed:
//			1) Specify Load Decrement
//			2) Specify Series Equivalent Load Parameters
//			3) Specify Parallel Equivalent Load Parameters
//			4) Specify Load Q
//		Comma Separated Variable (filename.csv) is also produced to allow easy
//		import into an MS Excel spreadsheet.
//
// Rev: j December 4, 2008
//
//		Allowed input of either Load Decrement or Q and Fraction Bandwidth.
//		If Q and Fraction Bandwidth is chosen, then Load Decrement is 
//		calculated from these values.
//		Load Decrement is utilized to obtain G(n) values in reference:
//			G. Matthaei, L. Young, and E.M.T. Jones,
//			Microwave Filters, Impedance Matching Networks, and
//			Coupling Structures,
//			Artech House, 1980.
//			See Figures 4.09-5 through 4.09-7 (pp 126-128).
//		This program offers a substitute to these figures.
//		Using Load Decrement as an input allows correlation of formulas
//		found in QEX article, "Broadband Impedance Matching", 
//		Frederic B. Hubler, QEX, November/December 2008, pp 23-29.
//		AL in BASIC program replaced by 1/LoadDecrement in this revision.
//		Refer to BASIC program listing at the end of this file.
//		LINE 210 AL = Q*BW ...
//		1/LoadDecrement = Q * (Fractional Bandwidth)
//		LoadDecrement = 1 / ( Q * (Fractional Bandwidth) )
//
// Rev: i November 26, 2008
//
//		Released for review by Frederic B. Hubler. 
//
//-----------------------------------------------------------------------------

#include <iostream>
#include <fstream>
#include <iomanip>
#include <math.h>
#include <stdio.h>
#include <time.h>

using namespace std;



//-----------------------------------------------------------------------------
//
// dPI  This returns pi.
double dPI (void)
{
	return ( acos( double(-1.0) ) );
}

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



//-----------------------------------------------------------------------------
//
// 100 DEF FNS
double FNS ( double X)
{
	// return ( (exp(X)-exp(-X))/2 );
	return ( sinh(X) );
}
//-----------------------------------------------------------------------------


//-----------------------------------------------------------------------------
//
// 110 DEF FNC
double FNC ( double X)
{
	// return ( (exp(X)+exp(-X))/2 );
	return ( cosh(X) );
}
//-----------------------------------------------------------------------------


//-----------------------------------------------------------------------------
//
// 120 DEF FNI
double FNI ( double X )
{
	return ( log(X+sqrt(X*X+double(1))) );
}
//-----------------------------------------------------------------------------


//-----------------------------------------------------------------------------
//
// 130 DEF FNG
double FNG ( double X, int N )
{
	return ( FNS(double(N)*X)/( FNC(double(N)*X) * FNC(X) ) );
}
//-----------------------------------------------------------------------------


//-----------------------------------------------------------------------------
//
// 140 DEF FNP
double FNP( double X, int N )
{
	return (	(double(N)-
				(  FNS(double(N)*X)*FNC(double(N)*X)*FNS(X)  ) / FNC(X) )
				/ 
				(  pow(FNC(double(N)*X),double(2))    * FNC(X) ) );
}
//-----------------------------------------------------------------------------


//-----------------------------------------------------------------------------
//
// 150 DEF FNV
double FNV ( double X )
{
	return ( (double(1)+X)/(double(1)-X) );
}
//-----------------------------------------------------------------------------


//-----------------------------------------------------------------------------
//
// 160 DEF FNN
double FNN ( int R, int N )
{
	return ( double(4) * sin( (R - double(0.5))*(dPI()/double(N)) ) * 
		     sin( (R + double(0.5))*(dPI()/double(N)) )
		   );
}
//-----------------------------------------------------------------------------


//-----------------------------------------------------------------------------
//
// 165 DEF FND
double FND ( double A, double B, int R, int N )
{
	double XX, YY, XY;

	XX = FNS ( A ) * FNS ( A );

	YY = FNS ( B ) * FNS ( B );

	XY = FNS ( A ) * FNS ( B );

	return ( XX + YY + 
		     pow(sin(double(R)*(dPI()/double(N))), double(2))-double(2) * 
			 XY *
			 (cos(double(R)*(dPI()/double(N))) )
		   );
}

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

//-----------------------------------------------------------------------------
//
// Reading a Whole Number 
// (Input string strW must have characters between '0' and '9' 
// with no other characters or else input will be requested again. )
//

int Read_Whole_Number ( void )
{

int			i,j;
int			testchar;

int			W = 3;
int			lenW;
double		dW, dWtot;
char		strW[40];
bool		Wgood = false;

	cin >> strW;
	lenW = strlen(strW);

	Wgood = true;
	for ( i = 0; i < lenW; i++ )
	{
		testchar = int(strW[i]);

		// character '0' = ASCII 48
		// character '9' = ASCII 57

		if ( !( ( testchar >= 48 ) && ( testchar <= 57 ) ) )
		{
			// character bad
			Wgood = false;
		}
	}


	if ( Wgood == true )
	{
		j = lenW - 1;
		dWtot = 0.0;
		for ( i = 0; i < lenW; i++ )
		{
			W = int(strW[i]) - 48;
			dW = double(W) * pow ( double (10.0), double (j) );
			dWtot = dWtot + dW;
			j--;
		}
	}
	else
	{
		dWtot = double(-1.0);
	}

	if ( dWtot < double(0.0) )
		return (int(-1));
	else
		return (int(dWtot));

}

//-----------------------------------------------------------------------------
//
// Reading a Fractional Number
// (Input string strF must have characters between '0' and '9' or '.' 
// but no other characters or else input will be requested again. )
//

double Read_Fractional_Number ( void )
{

int			i,j;
int			testchar;
int			lenF;

int			F;
double		dF, dFtot;
char		strF[40];
int			indexF_dec_pt;
bool		Fgood = false;

	cin >> strF;
	lenF = strlen(strF);

	Fgood = true;
	indexF_dec_pt = -1;  // indicates no decimal point found
	for ( i = 0; i < lenF; i++ )
	{
		testchar = int(strF[i]);

		// character '.' = ASCII 46
		// character '0' = ASCII 48
		// character '9' = ASCII 57

		if ( ( ( testchar >= 48 ) && ( testchar <= 57 ) ) )
		{
			// character good
			;
		}
		else if ( testchar == 46 )
		{
			if (indexF_dec_pt == -1)
			{
				// character good, first decimal point
				indexF_dec_pt = i;
			}
			else
			{
				// character bad, not first decimal point
				Fgood = false;
			}
		}
		else
		{
			// character bad
			Fgood = false;
		}
	}


	if ( Fgood == true )
	{
		if ( indexF_dec_pt == -1 )
		{
			// no decimal point	found					

			j = lenF - 1;
			dFtot = 0.0;
			for ( i = 0; i < lenF; i++ )
			{
				F = int(strF[i]) - 48;
				dF = double(F) * pow ( double(10), double(j) );
				dFtot = dFtot + dF;
				j--;
			}

		}
		else if ( indexF_dec_pt == 0 )
		{
			// decimal point in first character

			j = -1;
			dFtot = 0.0;
			for ( i = 1; i < lenF; i++ )
			{
				F = int(strF[i]) - 48;
				dF = double(F) * pow ( double(10), double(j) );
				dFtot = dFtot + dF;
				j--;
			}

		}

		else if ( indexF_dec_pt == (lenF-1) )
		{
			// decimal point in last character

			j = lenF - 2;
			dFtot = 0.0;
			for ( i = 0; i < (lenF-1); i++ )
			{
				F = int(strF[i]) - 48;
				dF = double(F) * pow( double(10), double(j) );
				dFtot = dFtot + dF;
				j--;
			}

		}

		else
		{
			// decimal point not in first nor in last character

			dFtot = 0.0;

			// whole decimal places

			j = indexF_dec_pt-1;

			for ( i = 0; i < indexF_dec_pt; i++ )
			{
				F = int(strF[i]) - 48;
				dF = double(F) * pow( double(10), double(j) );
				dFtot = dFtot + dF;
				j--;
			}

			// fractional decimal places

			j = -1;
			
			for ( i = indexF_dec_pt + 1; i < lenF; i++ )
			{
				F = int (strF[i]) - 48;
				dF = double(F) * pow( double(10), double(j) );
				dFtot = dFtot + dF;
				j--;
			}

		}
	}
	else
	{
		dFtot = double(-1);
	}

	if ( dFtot < double(0) )
		return (double(-1));
	else
		return (dFtot);

}

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

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

void main (void)
{


// begin code for main


// 160

const double dK1 =		double(1.7);
const double dK2 =		double(0.6);
const double dFVmin =	double(1.0E-9);

double		A, B;
double		ReflectCoefMag_Max, ReflectCoefMag_Min, SWR_Max, SWR_Min;
double		ReturnLoss_Max_dB, ReturnLoss_Min_dB;
double		MismatchLoss_Max_dB, MismatchLoss_Min_dB;
double		F1, F2, J1, J2, J3, J4, JD, DA, DB, FV;
int			IT, R;
int			i,j;

// 170
double		G[12], Gp[12], Gpp[12], L[12], C[12];
double		LoadResistance, SourceResistance;

// Input Variables
int			N;
double		LoadDecrement, Q, FractionalBandwidth;
char		inchar;
char		filename[40];
char		filename_txt[40];
char		filename_csv[40];
bool		continue_program = false;

double		deltaBW_ip, aveR_ip, deltaX_ip, aveG_ip, deltaB_ip;
double		Freq0, Freq1, Freq2, w0_Hz;
double		Xslope, Bslope;

enum		load_decrement_type_enum { load_decrement, load_parameters, load_q} load_decrement_type;
enum		serial_parallel_load_enum {serial, parallel} serial_parallel_load_type;

char		time_in[25];
char		time_out[25];

//-----------------------------------------------------------------------------
//
// Limits on Input Variables
//

//
// Reasonable limits for Load Decrement values based upon following: 
//
// Load Decrement values range from 0.1 to 2.0 in Figure 6.23 on page 198, in 
// text:
//		Circuit Design Using Personal Computers
//		Thomas R. Cuthbert, Jr.
//		1983 John Wiley and Sons 
//
// Load Decrement values range from 0.05 to 2.0 in Figures 4.09-5, 4.09-6, 
// 4.09-7, and 4.09-8 on page 198, in text:
//		G. Matthaei, L. Young, and E.M.T. Jones,
//		Microwave Filters, Impedance Matching Networks, and
//		Coupling Structures,
//		Artech House, 1980.
//

const double	LoadDecrement_min = double(0.01);
const double	LoadDecrement_max = double(20);

const double	Q_min = double(0.0333);
const double	Q_max = double(300);

const double	aveR_ip_min = double(0);
const double	aveR_ip_max = double(1000000);
const double	deltaX_ip_min = double(0);
const double	deltaX_ip_max = double(1000000);

const double	aveG_ip_min = double(0);
const double	aveG_ip_max = double(1000000);
const double	deltaB_ip_min = double(0);
const double	deltaB_ip_max = double(1000000);

const double	Freq_min = double(0.000001);	// MHz
const double	Freq_max = double(1000000);		// MHz

//-----------------------------------------------------------------------------
  
// Begin Major Loop
	
	do
	{

	//-----------------------------------------------------------------------------
	//
	// Get the current local time and date
	//

	time_t rawtime;
	struct tm * timeinfo;

	time ( &rawtime );
	timeinfo = localtime ( &rawtime );

	//-----------------------------------------------------------------------------
	//
	// Display to console the current local date and time
	//

	// cout	<< endl;
	// cout	<< "Current local date and time: ";
	// cout	<< asctime (timeinfo);

	strcpy(time_in,asctime (timeinfo));

	// Place four digit year after "Www Mmm dd" so output is:
	// "Www Mmm dd yyyy hh:mm:ss"
	for ( i = 0; i <= 10; i++)
		time_out[i] = time_in[i];

	for ( i = 11, j = 20; i <= 14; i++, j++)
		time_out[i] = time_in[j];

	for ( i = 15, j = 10; i <= 23; i++, j++)
		time_out[i] = time_in[j];

	time_out[24] = '\0';  // terminate string

	cout	<< endl;
	cout	<< "Current local date and time: " << time_out << endl;

	//-------------------------------------------------------------------------
	//
	// Header
	//

	cout	<< endl 
			<< "C++ version of Program B6-3" << endl
			<< "Newton - Raphson Solution then Levy Matching found in text:" << endl
			<< "Circuit Design Using Personal Computers" << endl
			<< "Thomas R. Cuthbert, Jr." << endl
			<< "1983 John Wiley and Sons" << endl
			<< endl
			<< "Cuthbert_B6_3n .cpp .exe" << endl
			<< "Rev: n January 30, 2009" << endl
			<< "Steven Schultz, WB8WGY" << endl
			<< endl;


	//-------------------------------------------------------------------------
	//
	// Input File Name
	//

	cout << "Enter Result File Name without suffix: ";
	cin >> filename;

	strcpy(filename_txt, filename);
	strcpy(filename_csv, filename);
	strcat(filename_txt, ".txt");
	strcat(filename_csv, ".csv");
	fstream fd_txt(filename_txt, ios::out);
	fstream fd_csv(filename_csv, ios::out);

	//-------------------------------------------------------------------------
	//
	// Input N
	//

	do
	{
		cout << endl
			 << "Enter N with a value >/= 1 and </= 10:  ";
		N = Read_Whole_Number ();
	} while ( !( (N >= 1) && (N <= 10 ) ) );

	//-------------------------------------------------------------------------
	//
	// Enter method to be used to obtain Load Decrement:
	//  1) Specify Load Decrement,
	//	2) Load Parameters,
	//  or
	//	3) Specify Load Q
	//

	do
	{

	cout << endl
		 << "Enter \"d\" or \"D\" to enter Load Decrement directly or" << endl
		 << "Enter \"l\" or \"L\" to obtain Load Decrement from Load Parameters or" << endl
		 << "Enter \"q\" or \"Q\" to obtain Load Decrement from Q and Fractional Bandwidth" << endl
		 << "Enter: ";
	cin  >> inchar;

	} while (  !(
				  (inchar == 'd') || (inchar == 'D') ||
				  (inchar == 'l') || (inchar == 'L') ||
				  (inchar == 'q') || (inchar == 'Q')
				)
			);

	switch (inchar)
	{

		case ('d'):
		case ('D'):
			load_decrement_type = load_decrement;
			break;

		case('l'):
		case('L'):
			load_decrement_type = load_parameters;
			break;

		case('q'):
		case('Q'):
			load_decrement_type = load_q;
			break;

	}

	//-------------------------------------------------------------------------
	//
	// Enter type of load:
	//	1) Series Load,
	//  or
	//	2) Parallel Load
	//

	do
	{

	cout << endl
		 << "Enter \"s\" or \"S\" for Series Load or Series Equivalent Load Parameters or" << endl
		 << "Enter \"p\" or \"P\" for Parallel Load or Parallel Equivalent Load Parameters" << endl
		 << "Enter: ";
	cin  >> inchar;

	} while (  !(
				  (inchar == 's') || (inchar == 'S') ||
				  (inchar == 'p') || (inchar == 'P')
				)
			);

	switch (inchar)
	{

		case('s'):
		case('S'):
			serial_parallel_load_type = serial;
			break;
		
		case ('p'):
		case ('P'):
			serial_parallel_load_type = parallel;
			break;

	}

	//-------------------------------------------------------------------------
	//
	// Enter remainder of inputs and
	// Begin Part 1: Preliminary
	//

	//-------------------------------------------------------------------------
	//
	// Enter f1 and f2.
	// Calculate f0, f2-f1, and wm.
	//

// 200

	cout << setprecision(6); 
	cout.setf(ios::fixed, ios::floatfield);

	do
	{
		cout << endl
			 << "Enter Lower Band Edge Frequency (f1) in MHz > " << Freq_min << " and < " << Freq_max << ":  ";  // 1Hz to 1000 GHz
		Freq1 = Read_Fractional_Number ();
	} while ( !( ( Freq1 > Freq_min ) && 
		         ( Freq1 < Freq_max )    ) );	

	do
	{
		cout << endl
			 << "Enter Higher Band Edge Frequency (f2) in MHz > " << Freq1 << " and < " << Freq_max << ":  ";  // Freq1 to 100 GHz
		Freq2 = Read_Fractional_Number ();
	} while ( !( ( Freq2 > Freq1    ) && 
			     ( Freq2 < Freq_max )    ) );

	Freq0 = sqrt( Freq1 * Freq2 );
	w0_Hz = double(2) * dPI() * Freq0 * pow ( double (10.0), double (6.0) ) ;

	FractionalBandwidth = ( Freq2 - Freq1 ) / Freq0; 

	cout << setprecision(6); 
	cout.setf(ios::fixed, ios::floatfield);

	cout	<< endl
			<< "Low Frequency        (f1)    (MHz): " << Freq1 << endl
			<< "Center Frequency     (f0)    (MHz): " << Freq0 << endl
			<< "High Frequency       (f2)    (MHz): " << Freq2 << endl
			<< "Match Bandwidth      (f2-f1) (MHz): " << (Freq2-Freq1) << endl
			<< "Fractional Match Bandwidth (wm)   : " << FractionalBandwidth << endl;

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

	if ( load_decrement_type != load_parameters )
	{

		if ( load_decrement_type == load_decrement )
		{
			
			cout << setprecision(3); 
			cout.setf(ios::fixed, ios::floatfield);

			// Input Load Decrement from console.
			do
			{
				cout << endl
					 << "Enter Load Decrement with a value >= " << LoadDecrement_min << " and <= " << LoadDecrement_max << ":  ";
				LoadDecrement = Read_Fractional_Number ();
			} while ( !( ( LoadDecrement >= LoadDecrement_min ) && 
						 ( LoadDecrement <= LoadDecrement_max )     ) 
					);
			
			Q = double(1) / ( LoadDecrement * FractionalBandwidth );
		}

		else if ( load_decrement_type == load_q )
		{

			cout << setprecision(5); 
			cout.setf(ios::fixed, ios::floatfield);

			cout	<< endl
					<< "Calculate Load Decrement from Q and Fractional Match Bandwidth ((f2-f1)/f0)" << endl;

			// Input Load Q from console.
			do
			{
			cout	<< endl
					<< "Enter Load Q with a value >= " << Q_min << " and <= " << Q_max << ":  ";
				Q = Read_Fractional_Number ();
			} while ( !( (Q >= Q_min) && (Q <= Q_max) ) );

			LoadDecrement = double(1) / (Q * FractionalBandwidth);
		}

		if ( serial_parallel_load_type == serial ) 
		{
			
			cout << setprecision(0); 
			cout.setf(ios::fixed, ios::floatfield);

			// Enter average resistance
			do
			{
			cout	<< endl
					<< "'average resistance' is used to impedance scale network values" << endl
					<< "'average resistance' is average series load resistance over the match band (f1 to f2)" << endl
					<< "Enter 'average resistance' with a value between " << aveR_ip_min << " and " << aveR_ip_max << " ohms: ";
				aveR_ip = Read_Fractional_Number ();
			} while ( !( ( aveR_ip >  aveR_ip_min ) && 
						 ( aveR_ip <= aveR_ip_max )    ) 
					);
		}
		else if ( serial_parallel_load_type == parallel )
		{
						
			cout << setprecision(0); 
			cout.setf(ios::fixed, ios::floatfield);

			// Enter average conductance
			do
			{
			cout	<< endl
					<< "'average conductance' is used to admittance scale network values" << endl
					<< "'average conductance' is average parallel load conductance over the match band (f1 to f2)" << endl
					<< "Enter 'average conductance' with a value between " << aveG_ip_min << " and " << aveG_ip_max << " Siemens: ";
				aveG_ip = Read_Fractional_Number ();
			} while ( !( ( aveG_ip >  aveG_ip_min ) && 
						 ( aveG_ip <= aveG_ip_max )    ) 
					);
		}

	}

	else if ( load_decrement_type == load_parameters )
	{

		if ( serial_parallel_load_type == serial ) 
		{
			// Series Load

			cout	<< endl
					<< "Calculate Load Decrement from the following series load slope parameters to be entered:"<< endl
					<< "   'average resistance' is average series load resistance over the match band (f1 to f2)" << endl
					<< "   'delta bandwidth' is bandwidth in MHz used to calculate 'delta reactance' " << endl
					<< "   'delta reactance' is change in load series reactance over 'delta bandwidth' evaluted at f0" << endl
					<< "   Note: 'delta bandwidth' is not 'match bandwidth'  " << endl;

			// Enter average resistance
			do
			{
			cout	<< endl
					<< "Enter 'average resistance' with a value between " << aveR_ip_min << " and " << aveR_ip_max << " ohms: ";
				aveR_ip = Read_Fractional_Number ();
			} while ( !( ( aveR_ip >  aveR_ip_min ) && 
						 ( aveR_ip <= aveR_ip_max )    ) );

			// Enter delta bandwidth
			do
			{
			cout	<< endl
					<< "Enter 'delta bandwidth' in MHz with a value between 0 and (" << (Freq2-Freq1) << "): ";
				deltaBW_ip = Read_Fractional_Number ();
			} while ( !( ( deltaBW_ip >  double(0)     ) && 
						 ( deltaBW_ip <= (Freq2-Freq1) )    )
					);

			// Enter delta reactance
			do
			{
			cout	<< endl
					<< "Enter 'delta reactance' with a value between " << deltaX_ip_min << " and "<< deltaX_ip_max << " ohms: ";
				deltaX_ip = Read_Fractional_Number ();
			} while ( !( ( deltaX_ip >  deltaX_ip_min ) && 
						 ( deltaX_ip <= deltaX_ip_max )     ) 
					);

			Xslope = ( Freq0 / double(2) ) * ( deltaX_ip / deltaBW_ip );
			LoadDecrement = ( aveR_ip / Xslope ) * ( 1 / FractionalBandwidth );

		}

		else if ( serial_parallel_load_type == parallel )
		{

			// Parallel Load

			cout	<< endl
					<< "Calculate Load Decrement from the following parallel load slope parameters to be entered:"<< endl
					<< "   'average conductance' is average parallel load conductance over the match band (f1 to f2)" << endl
					<< "   'delta bandwidth' is bandwidth in MHz used to calculate 'delta susceptance' " << endl
					<< "   'delta susceptance' is change in load parallel susceptance over 'delta bandwidth' evaluted at f0" << endl
					<< "   Note: 'delta bandwidth' is not 'match bandwidth' " << endl;

			// Enter average conductance
			do
			{
			cout	<< endl
					<< "Enter 'average conductance' with a value between " << aveG_ip_min << " and " << aveG_ip_max << " Siemens: ";
				aveG_ip = Read_Fractional_Number ();
			} while ( !( ( aveG_ip >  aveG_ip_min ) && 
						 ( aveG_ip <= aveG_ip_max )    ) 
					);

			// Enter delta bandwidth
			do
			{
			cout	<< endl
					<< "Enter 'delta bandwidth' in MHz with a value between 0 and (" << (Freq2-Freq1) << "): ";
				deltaBW_ip = Read_Fractional_Number ();
			} while ( !( ( deltaBW_ip >  double(0)     ) && 
						 ( deltaBW_ip <= (Freq2-Freq1) )    ) 
					);

			// Enter delta susceptance
			do
			{
			cout	<< endl
					<< "Enter 'delta susceptance' with a value between " << deltaB_ip_min << " and " << deltaB_ip_max << " Siemens: ";
				deltaB_ip = Read_Fractional_Number ();
			} while ( !( ( deltaB_ip >  deltaB_ip_min ) && 
						 ( deltaB_ip <= deltaB_ip_max )    )
					);

			Bslope = ( Freq0 / double(2) ) * ( deltaB_ip / deltaBW_ip );
			LoadDecrement = ( aveG_ip / Bslope ) * ( 1 / FractionalBandwidth );

		}
	
		Q = double(1) / ( LoadDecrement * FractionalBandwidth );

	}


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

	if ( LoadDecrement < LoadDecrement_min )
	{
		cout << endl;
		cout << "--------------------------------------------------------------------------------" << endl;
		cout << " WARNING: Load Decrement < Recommended minimun Load Decrement (" << LoadDecrement_min << ") " << endl;
		cout << "--------------------------------------------------------------------------------" << endl;
		cout << endl;
	}
	else if ( LoadDecrement > LoadDecrement_max )
	{
		cout << endl;
		cout << "--------------------------------------------------------------------------------" << endl;
		cout << " WARNING: Load Decrement > Recommended maximum Load Decrement (" << LoadDecrement_max << ") " << endl;
		cout << "--------------------------------------------------------------------------------" << endl;
		cout << endl;
	}

	if ( Q < Q_min )
	{
		cout << endl;
		cout << "--------------------------------------------------------------------------------" << endl;
		cout << " WARNING: Load Q < Recommended minimun Load Q (" << Q_min << ") " << endl;
		cout << "--------------------------------------------------------------------------------" << endl;
		cout << endl;
	}
	else if ( Q > Q_max )
	{
		cout << endl;
		cout << "--------------------------------------------------------------------------------" << endl;
		cout << " WARNING: Load Q > Recommended maximum Load Q (" << Q_max << ") " << endl;
		cout << "--------------------------------------------------------------------------------" << endl;
		cout << endl;
	}


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

	//
	// End Part 1
	//
	//-------------------------------------------------------------------------


	//-------------------------------------------------------------------------
	//
	// Begin Part 2: Fano's Broadband-Matching Limitations
	//				 Levy Matching
	//				 Newton - Raphson Solution
	//

	// 210
	// 300

	// 310
	// A = FNI(   (dK1*(pow(AL,dK2))+1) * (sin((dPI()/2)/N))  /AL );
	A = FNI(  (dK1*(pow(LoadDecrement,-dK2))+double(1)) * 
		      (sin((dPI()/double(2))/double(N))) * LoadDecrement);


	// 320
	// B = FNI(   (dK1*(pow(AL,dK2))-1) * (sin((dPI()/2)/N))  /AL );
    B = FNI(  (dK1*(pow(LoadDecrement,-dK2))-double(1)) * 
		      (sin((dPI()/double(2))/double(N))) * LoadDecrement);


	// 330
	if ( B < double(0) )
		B = double(0.001);

	// 350
	IT = 1;

	// 360

	// 400 REM: CALC FUNC, DELTA & NEW A,B
	do
	{

	// 410
	// F1=FNS(A)-FNS(B)-(2/AL)*sin((dPI()/2)/N);
	F1=FNS(A)-FNS(B)-(double(2)*LoadDecrement)*
		sin((dPI()/double(2))/double(N));

	// 420
	F2=FNG(A,N)-FNG(B,N);

	// 430
	J1=FNC(A);

	// 440
	J2=FNP(A,N);

	// 450
	J3=-FNC(B);

	// 460
	J4=-FNP(B,N);

	// 470
	JD=J1*J4-J2*J3;

	// 0.5 damping factor multiplied to DA and DB to ensure large values of Q and or BW
	// allow Newton - Raphson Solution to converge.

	// 480
	DA=(-(J4*F1-J3*F2)/JD)*double(0.5);	

	// 490
	DB=(-(-J2*F1+J1*F2)/JD)*double(0.5);

	// 500
	FV=F1*F1+F2*F2;

	// 510

		if ( FV > dFVmin )
		{
			// 530 IT=IT+1
			IT = IT + 1;

			// 540 A=A+DA
			A = A + DA;

			// 550 B=B+DB
			B = B + DB;
		}

	}
	// 520 IF FV<1.E-9 GOTO600
	// 560 GOTO400
	while ( FV > dFVmin );


// 600 PRINT"CONVERGERED"
	// Convergence criteria met (FV < dFVmin).

// 610 RH=FNC(N*B)/FNC(N*A)
	ReflectCoefMag_Max = FNC(double(N)*B)/FNC(double(N)*A);

// 620 RL=FNS(N*B)/FNS(N*A)
	ReflectCoefMag_Min = FNS(double(N)*B)/FNS(double(N)*A);

// 630 SH=FNV(RH)
	SWR_Max=FNV(ReflectCoefMag_Max);

// 640 SL=FNV(RL)
	SWR_Min=FNV(ReflectCoefMag_Min);

	ReturnLoss_Max_dB = double(20) * log10( ReflectCoefMag_Max );

	ReturnLoss_Min_dB = double(20) * log10( ReflectCoefMag_Min );

	MismatchLoss_Max_dB = double(10) * log10( 1 - pow( ReflectCoefMag_Max, double(2) ) );

	MismatchLoss_Min_dB = double(10) * log10( 1 - pow( ReflectCoefMag_Min, double(2) ) );

	//
	// End Part 2
	//
	//-------------------------------------------------------------------------


	//-------------------------------------------------------------------------
	//
	// Begin Part 3: Resistive Source with g(i) Prototype Values
	//

// 650 

// 800 REM:CALC ELEMENT VALUES

// 810 G(1)=AL
   // G[0] = AL;

// 820 FOR R=1 TO N-1
// 830 G(R+1)=FNN(R)/FND(R)/G(R)
// 840 NEXT R

// 850
// 860
// 870


// 874 GS=2/G(N)*SIN(PI/2/N)/(X+Y)
/*
(2*SIN(PI/2/N))
/
(G(N)*(X+Y))
*/

	// Calculate g(i) values

	G[0] = 1;
	G[1] = 1/LoadDecrement;

	for ( R = 1; R <= (N-1); R++)
	{
		G[R+1] = FNN( R, N ) / ( FND( A, B, R, N ) * G[R] );
	}

	G[N+1] = (double(2)*sin((dPI()/double(2))/double(N))) / 
		 (G[N]*(FNS(A)+FNS(B)));

	//
	// End Part 3
	//
	//-------------------------------------------------------------------------

// 860

	//-------------------------------------------------------------------------
	//
	// Begin Part 4: Network Element Values Calculated
	//

	if ( serial_parallel_load_type == serial ) 
	{
		//---------------------------------------------------------------------
		// Series Equivalent Load

		//---------------------------------------------------------------------
		// Scale to Load Resistance
		// ( aveR_ip )

		// shunt elements
		for ( R = 0; R <= (N+1); R=R+2)
		{
			Gp[R] = G[R] / aveR_ip;
		}

		// series elements
		for ( R = 1; R <= (N+1); R=R+2)
		{
			Gp[R] = G[R] * aveR_ip;
		}

		//---------------------------------------------------------------------
		// Scale to Frequency

		for ( R = 1; R <= (N); R=R+1)
		{
			Gpp[R] = Gp[R] / FractionalBandwidth;
		}

		//---------------------------------------------------------------------
		// Calculate L and C values from Center Frequency

		// series elements
		for ( R = 1; R <= (N); R=R+2)
		{
			L[R] =            Gpp[R]  / w0_Hz;
			C[R] = (double(1)/Gpp[R]) / w0_Hz;
		}

		// shunt elements
		for ( R = 2; R <= (N); R=R+2)
		{
			L[R] = (double(1)/Gpp[R]) / w0_Hz;
			C[R] =            Gpp[R]  / w0_Hz;
		}

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

		// Load Resistance
		// (shunt Load Resistance)
		LoadResistance = double(1)/ Gp[0];

		// Source Resistance
		if ( (N - ((N/2)*2)) == 1 )
			// N odd (shunt Source Resistance)
			SourceResistance = double(1)/ Gp[N+1];
		else
			// N even (series Source Resistance)
			SourceResistance = Gp[N+1];

	}
	else if ( serial_parallel_load_type == parallel )
	{
		//---------------------------------------------------------------------
		// Parallel Equivalent Load

		//---------------------------------------------------------------------
		// Scale to Load Resistance
		// ( 1 / aveG_ip )

		// series elements
		for ( R = 0; R <= (N+1); R=R+2)
		{
			Gp[R] = G[R] * ( double(1) / aveG_ip );
		}

		// shunt elements
		for ( R = 1; R <= (N+1); R=R+2)
		{
			Gp[R] = G[R] / ( double(1) / aveG_ip );
		}

		//---------------------------------------------------------------------
		// Scale to Frequency

		for ( R = 1; R <= (N); R=R+1)
		{
			Gpp[R] = Gp[R] / FractionalBandwidth;
		}

		//---------------------------------------------------------------------
		// Calculate L and C values from Center Frequency

		// shunt elements
		for ( R = 1; R <= (N); R=R+2)
		{
			L[R] = (double(1)/Gpp[R]) / w0_Hz;
			C[R] =            Gpp[R]  / w0_Hz;
		}

		// series elements
		for ( R = 2; R <= (N); R=R+2)
		{
			L[R] =            Gpp[R]  / w0_Hz;
			C[R] = (double(1)/Gpp[R]) / w0_Hz;
		}

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

		// Load Resistance
		// (series Load Resistance)
		LoadResistance = Gp[0];

		// Source Resistance
		if ( (N - ((N/2)*2)) == 1 )
			// N odd (series Source Resistance)
			SourceResistance = Gp[N+1];
		else
			// N even (shunt Source Resistance)
			SourceResistance = double(1) / Gp[N+1];
	}

	//
	// End Part 4
	//
	//-------------------------------------------------------------------------


//-----------------------------------------------------------------------------
//
//	Report Section
//
//-----------------------------------------------------------------------------


	//-----------------------------------------------------------------------------
	//
	//	Display to console
	//
	//-----------------------------------------------------------------------------

	//-------------------------------------------------------------------------
	//
	// Inputs
	//

	cout	<< endl;
	cout	<< "Inputs:" << endl;

	cout	<< "   Text Result File Name:                " << filename_txt  << endl;
	cout	<< "   CSV Result File Name:                 " << filename_csv  << endl;

	cout	<< setprecision(0); 
	cout.setf(ios::fixed, ios::floatfield);

	cout	<< "   N =                                   " << N << endl;

	cout	<< setprecision(6); 
	cout.setf(ios::fixed, ios::floatfield);

	cout	<< "   Low Frequency (f1) (MHz) =            " << Freq1 << endl
			<< "   High Frequency (f2) (MHz) =           " << Freq2 << endl;

	cout	<< setprecision(6);
	cout.setf(ios::scientific, ios::floatfield);

	switch (load_decrement_type)
	{

	case (load_decrement):
	cout	<< "   Load Decrement entered directly" << endl
			<< "      Load Decrement =                   " << LoadDecrement << endl;
			break;

	case (load_parameters):
	cout	<< "   Load Decrement obtained from Load Parameters" << endl;
			break;

	case (load_q):
	cout	<< "   Load Decrement obtained from Q and Fractional Bandwidth" << endl
			<< "      Load Q =                           " << Q << endl;
			break;

	}

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

	if ( load_decrement_type == load_parameters )
	{

	if (serial_parallel_load_type == serial)
	{
	cout	<< "      Series Load Parameters entered" << endl
			<< "         average resistance (ohms) =     " << aveR_ip << endl
			<< "         delta bandwidth (MHz) =         " << deltaBW_ip << endl
			<< "         delta reactance (ohms) =        " << deltaX_ip << endl;
	}
	else if (serial_parallel_load_type == parallel)
	{
	cout	<< "      Parallel Load Parameters entered" << endl
			<< "         average conductance (Siemens) = " << aveG_ip << endl
			<< "         delta bandwidth (MHz) =         " << deltaBW_ip << endl
			<< "         delta susceptance (Siemens) =   " << deltaB_ip << endl;
	}

	}
	else // load_decrement or load_q
	{

	if      (serial_parallel_load_type == serial)
	{
	cout	<< "      Series Load" << endl
			<< "         average resistance (ohms) =     " << aveR_ip << endl;
	}
	else if (serial_parallel_load_type == parallel)
	{
	cout	<< "      Parallel Load" << endl
			<< "         average conductance (Siemens) = " << aveG_ip << endl;
	}

	}

	//-------------------------------------------------------------------------
	//
	// Outputs
	//

	//-------------------------------------------------------------------------
	//
	// Output Part 1: Preliminary
	//

	cout	<< endl
			<< "Output Part 1: Preliminary" << endl;

	cout	<< setprecision(6); 
	cout.setf(ios::fixed, ios::floatfield);

	cout	<< "   Center Frequency     (f0)    (MHz) =      " << Freq0 << endl
			<< "   Match Bandwidth      (f2-f1) (MHz) =      " << (Freq2-Freq1) << endl;

	cout	<< setprecision(9); 
	cout.setf(ios::fixed, ios::floatfield);

	cout	<< "   Fractional Match Bandwidth ((f2-f1)/f0) (wm) =  " << FractionalBandwidth << endl;

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

	cout << setprecision(9);
	cout.setf(ios::scientific, ios::floatfield);

	if		( load_decrement_type == load_parameters )
	{

	if (serial_parallel_load_type == serial)
	{
	cout	<< "   Reactance Slope Parameter =               " << Xslope << endl;
	}
	else if (serial_parallel_load_type == parallel)
	{
	cout	<< "   Susceptance Slope Parameter =             " << Bslope << endl;
	}
	
	}
	else if ( load_decrement_type == load_decrement )
	{

	if      (serial_parallel_load_type == serial)
	{
	cout	<< "   No Reactance Slope Parameter Calculated because Load Decrement entered" << endl;
	}
	else if (serial_parallel_load_type == parallel)
	{
	cout	<< "   No Susceptance Slope Parameter Calculated because Load Decrement entered" << endl;
	}

	}
	else if ( load_decrement_type == load_q )
	{

	if      (serial_parallel_load_type == serial)
	{
	cout	<< "   No Reactance Slope Parameter Calculated because Load Q entered" << endl;
	}
	else if (serial_parallel_load_type == parallel)
	{
	cout	<< "   No Susceptance Slope Parameter Calculated because Load Q entered" << endl;
	}

	}

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

	switch (load_decrement_type)
	{

	case (load_decrement):
	cout	<< "   Load Q =                                  " << Q << endl;
			break;

	case (load_parameters):
	cout	<< "   Load Decrement =                          " << LoadDecrement << endl;
	cout	<< "   Load Q =                                  " << Q << endl;
			break;

	case (load_q):
	cout	<< "   Load Decrement =                          " << LoadDecrement << endl;
			break;

	}

	//-------------------------------------------------------------------------
	//
	// Output Part 2: Fano's Broadband-Matching Limitations
	//				  Levy Matching
	//				  Newton - Raphson Solution
	//

	cout	<< endl;
	cout	<< "Output Part 2: Fano; Levy; Newton - Raphson Solution" << endl;

	cout << setprecision(3);
	cout.setf(ios::scientific, ios::floatfield);

	cout	<< "   FV =                   " << FV << endl;

	cout << setprecision(0);
	cout.setf(ios::fixed, ios::floatfield);

	cout	<< "   IT (iterations in Newton - Raphson Solution) = " << IT << endl;

	cout << setprecision(9);
	cout.setf(ios::fixed, ios::floatfield);

	cout 	<< "   A  =                   " << A << endl
			<< "   B  =                   " << B << endl
			<< "   Magnitude of Reflection Coefficient Min = " << ReflectCoefMag_Min << endl
			<< "   Magnitude of Reflection Coefficient Max = " << ReflectCoefMag_Max << endl
			<< "   Mismatch Loss Min (dB) = " << MismatchLoss_Min_dB << endl
			<< "   Mismatch Loss Max (dB) = " << MismatchLoss_Max_dB << endl
			<< "   Return Loss Min (dB) = " << ReturnLoss_Min_dB << endl
			<< "   Return Loss Max (dB) = " << ReturnLoss_Max_dB << endl
			<< "   SWR Min =              " << SWR_Min << endl
			<< "   SWR Max =              " << SWR_Max << endl;

	//-------------------------------------------------------------------------
	//
	// Output Part 3: Resistive Source with g(i) Prototype Values
	//

	cout	<< endl;
	cout	<< "Output Part 3: Resistive Source with g(i) Prototype Values" << endl;

	R = 0;
	cout	<< "   G(" << R << ") (load)   =     " << G[R] << endl;

	for ( R = 1; R <= N; R++)
	{
	cout	<< "   G(" << R << ")          =     " << G[R] << endl;
	}

	R = (N+1);
	cout	<< "   G(" << R << ") (source) =     " << G[R] << endl;

	//-------------------------------------------------------------------------
	//
	// Output Part 4: Network Element Values Calculated
	//

	cout	<< endl;
	cout	<< "Output Part 4: Network Element Values Calculated" << endl;

	cout << setprecision(9);
	cout.setf(ios::scientific, ios::floatfield);

	cout	<< "   Load Resistance (ohms) =   " << LoadResistance << endl;

	for ( R = 1; R <= N; R++)
	{
	cout	<< "   L(" << R << ") (H) =                 " << L[R] << endl;
	cout	<< "   C(" << R << ") (F) =                 " << C[R] << endl;
	}

	cout	<< "   Source Resistance (ohms) = " << SourceResistance << endl;

	//-------------------------------------------------------------------------
	//
	// Write to Text Result File
	//
	//-------------------------------------------------------------------------

	//-----------------------------------------------------------------------------
	//
	// Current local time and date
	//

	// fd_txt	<< endl;
	// fd_txt	<< "Current local time and date: ";
	// fd_txt	<< asctime (timeinfo);

	fd_txt	<< endl;
	fd_txt	<< "Current local date and time: " << time_out << endl;

	//-------------------------------------------------------------------------
	//
	// Header
	//

	fd_txt	<< endl 
			<< "C++ version of Program B6-3" << endl
			<< "Newton - Raphson Solution then Levy Matching found in text:" << endl
			<< "Circuit Design Using Personal Computers" << endl
			<< "Thomas R. Cuthbert, Jr." << endl
			<< "1983 John Wiley and Sons" << endl
			<< endl
			<< "Cuthbert_B6_3n .cpp .exe" << endl
			<< "Rev: n January 30, 2009" << endl
			<< "Steven Schultz, WB8WGY" << endl
			<< endl;

	fd_txt	<< setprecision(6); 
	fd_txt.setf(ios::fixed, ios::floatfield);

	//-------------------------------------------------------------------------
	//
	// Inputs
	//

	fd_txt	<< "Inputs:" << endl;

	fd_txt	<< "   Text Result File Name:                " << filename_txt << endl;
	fd_txt	<< "   CSV Result File Name:                 " << filename_csv << endl;

	fd_txt	<< setprecision(0); 
	fd_txt.setf(ios::fixed, ios::floatfield);

	fd_txt	<< "   N =                                   " << N << endl;

	fd_txt	<< setprecision(6); 
	fd_txt.setf(ios::fixed, ios::floatfield);

	fd_txt	<< "   Low Frequency (f1) (MHz) =            " << Freq1 << endl
			<< "   High Frequency (f2) (MHz) =           " << Freq2 << endl;

	fd_txt	<< setprecision(6);
	fd_txt.setf(ios::scientific, ios::floatfield);

	switch (load_decrement_type)
	{

	case (load_decrement):
	fd_txt	<< "   Load Decrement entered directly" << endl
			<< "      Load Decrement =                   " << LoadDecrement << endl;
			break;

	case (load_parameters):
	fd_txt	<< "   Load Decrement obtained from Load Parameters" << endl;
			break;

	case (load_q):
	fd_txt	<< "   Load Decrement obtained from Q and Fractional Bandwidth" << endl
			<< "      Load Q =                           " << Q << endl;
			break;

	}

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

	if ( load_decrement_type == load_parameters )
	{

	if (serial_parallel_load_type == serial)
	{
	fd_txt	<< "      Series Load Parameters entered" << endl
			<< "         average resistance (ohms) =     " << aveR_ip << endl
			<< "         delta bandwidth (MHz) =         " << deltaBW_ip << endl
			<< "         delta reactance (ohms) =        " << deltaX_ip << endl;
	}
	else if (serial_parallel_load_type == parallel)
	{
	fd_txt	<< "      Parallel Load Parameters entered" << endl
			<< "         average conductance (Siemens) = " << aveG_ip << endl
			<< "         delta bandwidth (MHz) =         " << deltaBW_ip << endl
			<< "         delta susceptance (Siemens) =   " << deltaB_ip << endl;
	}

	}
	else // load_decrement or load_q
	{

	if      (serial_parallel_load_type == serial)
	{
	fd_txt	<< "      Series Load" << endl
			<< "         average resistance (ohms) =     " << aveR_ip << endl;
	}
	else if (serial_parallel_load_type == parallel)
	{
	fd_txt	<< "      Parallel Load" << endl
			<< "         average conductance (Siemens) = " << aveG_ip << endl;
	}

	}

	//-------------------------------------------------------------------------
	//
	// Outputs
	//

	//-------------------------------------------------------------------------
	//
	// Output Part 1: Preliminary
	//

	fd_txt	<< endl
			<< "Output Part 1: Preliminary" << endl;

	fd_txt	<< setprecision(6); 
	fd_txt.setf(ios::fixed, ios::floatfield);

	fd_txt	<< "   Center Frequency     (f0)    (MHz) =      " << Freq0 << endl
			<< "   Match Bandwidth      (f2-f1) (MHz) =      " << (Freq2-Freq1) << endl;

	fd_txt	<< setprecision(9); 
	fd_txt.setf(ios::fixed, ios::floatfield);

	fd_txt	<< "   Fractional Match Bandwidth ((f2-f1)/f0) (wm) =  " << FractionalBandwidth << endl;

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

	fd_txt	<< setprecision(9);
	fd_txt.setf(ios::scientific, ios::floatfield);

	if		( load_decrement_type == load_parameters )
	{

	if (serial_parallel_load_type == serial)
	{
	fd_txt	<< "   Reactance Slope Parameter =               " << Xslope << endl;
	}
	else if (serial_parallel_load_type == parallel)
	{
	fd_txt	<< "   Susceptance Slope Parameter =             " << Bslope << endl;
	}
	
	}
	else if ( load_decrement_type == load_decrement )
	{

	if      (serial_parallel_load_type == serial)
	{
	fd_txt	<< "   No Reactance Slope Parameter Calculated because Load Decrement entered" << endl;
	}
	else if (serial_parallel_load_type == parallel)
	{
	fd_txt	<< "   No Susceptance Slope Parameter Calculated because Load Decrement entered" << endl;
	}

	}
	else if ( load_decrement_type == load_q )
	{

	if      (serial_parallel_load_type == serial)
	{
	fd_txt	<< "   No Reactance Slope Parameter Calculated because Load Q entered" << endl;
	}
	else if (serial_parallel_load_type == parallel)
	{
	fd_txt	<< "   No Susceptance Slope Parameter Calculated because Load Q entered" << endl;
	}

	}

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

	switch (load_decrement_type)
	{

	case (load_decrement):
	fd_txt	<< "   Load Q =                                  " << Q << endl;
			break;

	case (load_parameters):
	fd_txt	<< "   Load Decrement =                          " << LoadDecrement << endl;
	fd_txt	<< "   Load Q =                                  " << Q << endl;
			break;

	case (load_q):
	fd_txt	<< "   Load Decrement =                          " << LoadDecrement << endl;
			break;

	}

	//-------------------------------------------------------------------------
	//
	// Output Part 2: Fano's Broadband-Matching Limitations
	//				  Levy Matching
	//				  Newton - Raphson Solution
	//

	fd_txt	<< endl;
	fd_txt	<< "Output Part 2: Fano; Levy; Newton - Raphson Solution" << endl;

	fd_txt	<< setprecision(3);
	fd_txt.setf(ios::scientific, ios::floatfield);

	fd_txt	<< "   FV =                   " << FV << endl;

	fd_txt	<< setprecision(0);
	fd_txt.setf(ios::fixed, ios::floatfield);

	fd_txt	<< "   IT (iterations in Newton - Raphson Solution) = " << IT << endl;

	fd_txt	<< setprecision(9);
	fd_txt.setf(ios::fixed, ios::floatfield);

	fd_txt 	<< "   A  =                   " << A << endl
			<< "   B  =                   " << B << endl
			<< "   Magnitude of Reflection Coefficient Min = " << ReflectCoefMag_Min << endl
			<< "   Magnitude of Reflection Coefficient Max = " << ReflectCoefMag_Max << endl
			<< "   Mismatch Loss Min (dB) = " << MismatchLoss_Min_dB << endl
			<< "   Mismatch Loss Max (dB) = " << MismatchLoss_Max_dB << endl
			<< "   Return Loss Min (dB) = " << ReturnLoss_Min_dB << endl
			<< "   Return Loss Max (dB) = " << ReturnLoss_Max_dB << endl
			<< "   SWR Min =              " << SWR_Min << endl
			<< "   SWR Max =              " << SWR_Max << endl;

	//-------------------------------------------------------------------------
	//
	// Output Part 3: Resistive Source with g(i) Prototype Values
	//

	fd_txt	<< endl;
	fd_txt	<< "Output Part 3: Resistive Source with g(i) Prototype Values" << endl;

	R = 0;
	fd_txt	<< "   G(" << R << ") (load)   =     " << G[R] << endl;

	for ( R = 1; R <= N; R++)
	{
	fd_txt	<< "   G(" << R << ")          =     " << G[R] << endl;
	}

	R = (N+1);
	fd_txt	<< "   G(" << R << ") (source) =     " << G[R] << endl;

	//-------------------------------------------------------------------------
	//
	// Output Part 4: Network Element Values Calculated
	//

	fd_txt	<< endl;
	fd_txt	<< "Output Part 4: Network Element Values Calculated" << endl;

	fd_txt	<< setprecision(9);
	fd_txt.setf(ios::scientific, ios::floatfield);

	fd_txt	<< "   Load Resistance (ohms) =   " << LoadResistance << endl;

	for ( R = 1; R <= N; R++)
	{
	fd_txt	<< "   L(" << R << ") (H) =                 " << L[R] << endl;
	fd_txt	<< "   C(" << R << ") (F) =                 " << C[R] << endl;
	}

	fd_txt	<< "   Source Resistance (ohms) = " << SourceResistance << endl;


	//-----------------------------------------------------------------------------
	//
	//	Write to comma separated variable file
	//
	//-----------------------------------------------------------------------------

	//-----------------------------------------------------------------------------
	//
	// Current local time and date
	//

	// fd_csv	<< endl;
	// fd_csv	<< "Current local time and date: ";
	// fd_csv	<< asctime (timeinfo);

	fd_csv	<< endl;
	fd_csv	<< "Current local date and time: " << time_out << endl;

	//-------------------------------------------------------------------------
	//
	// Header
	//

	fd_csv	<< endl
			<< "C++ version of Program B6-3" << endl
			<< "Newton - Raphson Solution then Levy Matching found in text:" << endl
			<< "Circuit Design Using Personal Computers" << endl
			<< "Thomas R. Cuthbert Jr." << endl
			<< "1983 John Wiley and Sons" << endl
			<< endl
			<< "Cuthbert_B6_3n .cpp .exe" << endl
			<< "Rev: n January 30 2009" << endl
			<< "Steven Schultz WB8WGY" << endl
			<< endl;

	fd_csv	<< setprecision(6); 
	fd_csv.setf(ios::fixed, ios::floatfield);

	//-------------------------------------------------------------------------
	//
	// Inputs
	//

	fd_csv	<< "Inputs:" << endl;

	fd_csv	<< "   Text Result File Name:                ," << filename_txt << endl;
	fd_csv	<< "   CSV Result File Name:                 ," << filename_csv << endl;

	fd_csv	<< setprecision(0); 
	fd_csv.setf(ios::fixed, ios::floatfield);

	fd_csv	<< "   N =                                   ," << N << endl;

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

	switch (load_decrement_type)
	{
		
		case (load_decrement):
			fd_csv	<< "   D ,Load Decrement entered directly" << endl;
			break;
		
		case(load_parameters):
			fd_csv	<< "   L ,Load Decrement from Load Parameters" << endl;
			break;

		case(load_q):
			fd_csv	<< "   Q ,Load Decrement from Q and Fractional Bandwidth" << endl;
			break;

	}

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

	switch (serial_parallel_load_type)
	{
		case(serial):
			fd_csv	<< "   S ,Series Load or Series Equivalent Load Parameters" << endl;
			break;

		case (parallel):
			fd_csv	<< "   P ,Parallel Load or Parallel Equivalent Load Parameters" << endl;
			break;
	}

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

	fd_csv	<< setprecision(6); 
	fd_csv.setf(ios::fixed, ios::floatfield);

	fd_csv	<< "   Low Frequency (f1) (MHz) =            ," << Freq1 << endl
			<< "   High Frequency (f2) (MHz) =           ," << Freq2 << endl;

	fd_csv	<< setprecision(6);
	fd_csv.setf(ios::scientific, ios::floatfield);

	switch (load_decrement_type)
	{

	case (load_decrement):
	fd_csv	<< "   Load Decrement entered directly" << endl
			<< "      Load Decrement =                   ," << LoadDecrement << endl
			<< "                                               " << endl;
			break;

	case (load_parameters):
	fd_csv	<< "   Load Decrement obtained from Load Parameters" << endl
			<< "                                               " << endl
			<< "                                               " << endl;
			break;

	case (load_q):
	fd_csv	<< "   Load Decrement obtained from Q and Fractional Bandwidth" << endl
			<< "                                               " << endl
			<< "      Load Q =                           ," << Q << endl;
			break;

	}

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

	if ( load_decrement_type == load_parameters )
	{

	if (serial_parallel_load_type == serial)
	{
	fd_csv	<< "      Series Load Parameters entered" << endl
			<< "         average resistance (ohms) =     ," << aveR_ip << endl
			<< "         delta bandwidth (MHz) =         ," << deltaBW_ip << endl
			<< "         delta reactance (ohms) =        ," << deltaX_ip << endl;
	}
	else if (serial_parallel_load_type == parallel)
	{
	fd_csv	<< "      Parallel Load Parameters entered" << endl
			<< "         average conductance (Siemens) = ," << aveG_ip << endl
			<< "         delta bandwidth (MHz) =         ," << deltaBW_ip << endl
			<< "         delta susceptance (Siemens) =   ," << deltaB_ip << endl;
	}

	}
	else // load_decrement or load_q
	{

	if      (serial_parallel_load_type == serial)
	{
	fd_csv	<< "      Series Load" << endl
			<< "         average resistance (ohms) =     ," << aveR_ip << endl
			<< "                                         " << endl
			<< "                                         " << endl;
	}
	else if (serial_parallel_load_type == parallel)
	{
	fd_csv	<< "      Parallel Load" << endl
			<< "         average conductance (Siemens) = ," << aveG_ip << endl
			<< "                                         " << endl
			<< "                                         " << endl;
	}

	}

	//-------------------------------------------------------------------------
	//
	// Outputs
	//

	//-------------------------------------------------------------------------
	//
	// Output Part 1: Preliminary
	//

	fd_csv	<< endl
			<< "Output Part 1: Preliminary" << endl;

	fd_csv	<< setprecision(6); 
	fd_csv.setf(ios::fixed, ios::floatfield);

	fd_csv	<< "   Center Frequency     (f0)    (MHz) =      ," << Freq0 << endl
			<< "   Match Bandwidth      (f2-f1) (MHz) =      ," << (Freq2-Freq1) << endl;

	fd_csv	<< setprecision(9); 
	fd_csv.setf(ios::fixed, ios::floatfield);

	fd_csv	<< "   Fractional Match Bandwidth ((f2-f1)/f0) (wm) =  ," << FractionalBandwidth << endl;

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

	fd_csv	<< setprecision(9);
	fd_csv.setf(ios::scientific, ios::floatfield);

	if		( load_decrement_type == load_parameters )
	{

	if (serial_parallel_load_type == serial)
	{
	fd_csv	<< "   Reactance Slope Parameter =               ," << Xslope << endl;
	}
	else if (serial_parallel_load_type == parallel)
	{
	fd_csv	<< "   Susceptance Slope Parameter =             ," << Bslope << endl;
	}
	
	}
	else if ( load_decrement_type == load_decrement )
	{

	if      (serial_parallel_load_type == serial)
	{
	fd_csv	<< "   No Reactance Slope Parameter Calculated because Load Decrement entered" << endl;
	}
	else if (serial_parallel_load_type == parallel)
	{
	fd_csv	<< "   No Susceptance Slope Parameter Calculated because Load Decrement entered" << endl;
	}

	}
	else if ( load_decrement_type == load_q )
	{

	if      (serial_parallel_load_type == serial)
	{
	fd_csv	<< "   No Reactance Slope Parameter Calculated because Load Q entered" << endl;
	}
	else if (serial_parallel_load_type == parallel)
	{
	fd_csv	<< "   No Susceptance Slope Parameter Calculated because Load Q entered" << endl;
	}

	}

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

	switch (load_decrement_type)
	{

	case (load_decrement):
	fd_csv	<< "   Load Decrement entered at input           " << endl
			<< "   Load Q =                                 ," << Q << endl;
			break;

	case (load_parameters):
	fd_csv	<< "   Load Decrement =                         ," << LoadDecrement << endl
			<< "   Load Q =                                 ," << Q << endl;
			break;

	case (load_q):
	fd_csv	<< "   Load Decrement =                         ," << LoadDecrement << endl
			<< "   Load Q entered at input                   " << endl;
			break;

	}

	//-------------------------------------------------------------------------
	//
	// Output Part 2: Fano's Broadband-Matching Limitations
	//				  Levy Matching
	//				  Newton - Raphson Solution
	//

	fd_csv	<< endl;
	fd_csv	<< "Output Part 2: Fano; Levy; Newton - Raphson Solution" << endl;

	fd_csv	<< setprecision(3);
	fd_csv.setf(ios::scientific, ios::floatfield);

	fd_csv	<< "   FV =                   ," << FV << endl;

	fd_csv	<< setprecision(0);
	fd_csv.setf(ios::fixed, ios::floatfield);

	fd_csv	<< "   IT (iterations in Newton - Raphson Solution) = ," << IT << endl;

	fd_csv	<< setprecision(9);
	fd_csv.setf(ios::fixed, ios::floatfield);

	fd_csv 	<< "   A  =                   ," << A << endl
			<< "   B  =                   ," << B << endl
			<< "   Magnitude of Reflection Coefficient Min = ," << ReflectCoefMag_Min << endl
			<< "   Magnitude of Reflection Coefficient Max = ," << ReflectCoefMag_Max << endl
			<< "   Mismatch Loss Min (dB) = ," << MismatchLoss_Min_dB << endl
			<< "   Mismatch Loss Max (dB) = ," << MismatchLoss_Max_dB << endl
			<< "   Return Loss Min (dB) = ," << ReturnLoss_Min_dB << endl
			<< "   Return Loss Max (dB) = ," << ReturnLoss_Max_dB << endl
			<< "   SWR Min =              ," << SWR_Min << endl
			<< "   SWR Max =              ," << SWR_Max << endl;

	//-------------------------------------------------------------------------
	//
	// Output Part 3: Resistive Source with g(i) Prototype Values
	//

	fd_csv	<< endl;
	fd_csv	<< "Output Part 3: Resistive Source with g(i) Prototype Values" << endl;

	R = 0;
	fd_csv	<< "   G(" << R << ") (load)   =     ," << G[R] << endl;

	for ( R = 1; R <= N; R++)
	{
	fd_csv	<< "   G(" << R << ")          =     ," << G[R] << endl;
	}

	R = (N+1);
	fd_csv	<< "   G(" << R << ") (source) =     ," << G[R] << endl;

	//-------------------------------------------------------------------------
	//
	// Output Part 4: Network Element Values Calculated
	//

	fd_csv	<< endl;
	fd_csv	<< "Output Part 4: Network Element Values Calculated" << endl;

	fd_csv	<< setprecision(9);
	fd_csv.setf(ios::scientific, ios::floatfield);

	fd_csv	<< "   Load Resistance (ohms) =   ," << LoadResistance << endl;

	for ( R = 1; R <= N; R++)
	{
	fd_csv	<< "   L(" << R << ") (H) =                 ," << L[R] << endl;
	fd_csv	<< "   C(" << R << ") (F) =                 ," << C[R] << endl;
	}

	fd_csv	<< "   Source Resistance (ohms) = ," << SourceResistance << endl;


	//-----------------------------------------------------------------------------
	//
	// Close files and determine if another pass is requested.
	//
	//-----------------------------------------------------------------------------

	fd_txt.close();
	fd_csv.close();

		do
		{

		cout << endl
			 << "Enter \"e\" or \"E\" to exit  " << endl
			 << "Enter \"c\" or \"C\" to continue  " << endl;
		cout << "Enter: ";
		cin  >> inchar;

		} while (  !(
						(inchar == 'e') || (inchar == 'E') ||
						(inchar == 'c') || (inchar == 'C')
				    )
				);


		if ( (inchar == 'c') || (inchar == 'C') )
			continue_program = true;
		else if ( (inchar == 'e') || (inchar == 'E') )
			continue_program = false;
		else
			continue_program = false;

	} while (continue_program);

	return;


} // end of main

// do not re-enter program
// 880 GOTO200
// 890 END


//*****************************************************************************
//
//
// Scanned from book:
// "Circuit Design Using Personal Computers"
// by Thomas R. Cuthbert, Jr.
// 1983 John Wiley and Sons, Inc
// 
// Scanned and Checked Oct. 30, 2008
//
// BASIC Language Program B6-3, Appendix B, page 443
// See section 6.4 on page 200 for explanation.
// Program B6-3. Levy Matching to Resistive Source with gi Prototype Values
// 10 	REM LEVY RS WITH G(N+1) MODIFICATION.
// 20 	REM 810116. TRC.
// 100 	DEF FNS(X)=(EXP(X)-EXP(-X))/2
// 110 	DEF FNC(X)=(EXP(X)+EXP(-X))/2 
// 120 	DEF FNI(X)=LOG(X+SQR(X*X+1))
// 130 	DEF FNG(X)=FNS(N*X)/FNC(N*X)/FNC(X) 
// 140 	DEF FNP(X)=(N-(FNS(N*X)*FNC(N*X)*FNS(X))/FNC(X))/(FNC(N*X)^2)/FNC(X)
// 150 	DEF FNV(X)=(1+X)/(1-X) 
// 160 	DEF FNN(R)=4*SIN((R-.5)*PN)*SIN((R+.5)*PN) 
// 165 	DEF FND(R)=X*X+Y*Y+SIN(R*PN)^2-2*X*Y*COS(R*PN)
// 170 	DIM G(10) 
// 180 	PI=3.1415926
// 200 	PRINT 
// 205 	PRINT"N,QL,%BW=";:INPUT N,Q,SW
// 210 	AL=Q*BW/l00
// 300 	REM: CALC INIT A,B BY CUTHBERT - 
// 310 	A=FNI((1.7*AL^.6+1)*(SIN(PI/2/N))/AL)
// 320 	B=FNI((1.7*AL^.6-1)*(SIN(PI/2/N))/AL)
// 330 	IF B<O THEN B=.OOl
// 350 	IT=O 
// 360 	PN=PI/N
// 400 	REM: CALC FUNC, DELTA & NEW A,B
// 410 	F1=FNS(A)-FNS(B)-(2/AL)*SIN(PI/2/N)
// 420 	F2=FNG(A)-FNG(B)
// 430 	J1=FNC(A)
// 440 	J2=FNP(A)
// 450 	J3=-FNC(B)
// 460 	J4=-FNP(B)
// 470 	JD=J1*J4-J2*J3
// 480 	DA=-(J4*F1-J3*F2)/JD 
// 490 	DB=-(-J2*F1+J1*F2)/JD 
// 500 	FV=F1*F1+F2*F2
// 520 	IF FV<1.E-9 GOT0600
// 530 	IT=IT+1
// 540 	A=A+DA
// 550 	B=B+DB
// 560 	GOT0400
// 600 	REM:CALC REFLECTION & SWR EXTREMES
// 602 	X=FNS(A)
// 604 	Y=FNS(B)
// 610 	RH=FNC(N*B)/FNC(N*A)
// 620 	RL=FNS(N*B)/FNS(N*A)
// 630 	SH=FNV(RH)
// 640 	SL=FNV(RL)
// 650 	PRINT"SWR FROM";SL;"TO";SH
// 720 	IF N>10 GOT0200
// 800 	REM:CALC ELEMENT VALUES
// 810 	G(1)=AL
// 820 	FOR R=1 TO N-1
// 830 	G(R+1)=FNN(R)/FND(R)/G(R)
// 840 	NEXT R
// 850 	FOR R=1 TO N
// 860 	PRINT"G(";R;")=";G(R)
// 870 	NEXT R
// 874 	GS=2/G(N)*SIN(PI/2/N)/(X+Y)
// 876 	PRINT"G(";N+l;")=";GS
// 880 	GOT0200
// 890 	END
