/* 
 * METRIC-Application:
 * four layer slab waveguide, 
 * leaky mode analysis 
 */

#include<stdio.h>
#include<stdlib.h>
#include<math.h>
#include"metric.h"

Complex Wgpnc(1.00, 0.0); // cover: air
Complex Wgpnf(3.45, 0.0); // film refractive index
Complex Wgpnb(1.45, 0.0); // buffer layer, refractive index
Complex Wgpns(3.45, 0.0); // substrate refractive index
double  Wgpt  =  0.22;    // slab thickness 
double  Wgpb  =  0.5;     // buffer thickness 
double  Wavel =  1.55;    // vacuum wavelength 

/* waveguide definition */
cWaveguide wgdef()
{
	// a waveguide with two inner layers
	cWaveguide g(2);

	// positions of the dielectric interfaces
	g.hx(2) = Wgpt;
	g.hx(1) = 0.0;
	g.hx(0) = -Wgpb;

	// refractive index values 
	g.n(3) = Wgpnc;
	g.n(2) = Wgpnf;
	g.n(1) = Wgpnb;
	g.n(0) = Wgpns;

	// the vacuum wavelength
	g.lambda = Wavel;

	return g;
}

/* calculate the leaky modes, write profile data */
int main()
{
	cModeArray ma;
		
	// define the waveguide
	cWaveguide wg = wgdef();

	// display window
	Interval disp(wg.hx(0)-2.0, wg.hx(wg.nx)+2.0);
	// inspect the refractive index profile
	wg.plot(disp, '0', '0');

	// find leaky TE modes
	leakymodeanalysis(wg, TE, ma);

	// mode data & profile plots, principal component 
	for(int i=0; i<=ma.num-1; ++i)
	{
		fprintf(stderr, "[%d] %s: neff = ", i, ma(i).ids); ma(i).neff.nice(stderr); fprintf(stderr, ", Lp = %g\n", ma(i).Lp()); 
		ma(i).plot(EY, REP, disp, 500, dig10(i), dig1(i), 'L');
		ma(i).plot(EY, IMP, disp, 500, dig10(i), dig1(i), 'V');
		ma(i).plot(EY, MOD, disp, 500, dig10(i), dig1(i), 'L');
	}

	// ... the same for the TM fields
	leakymodeanalysis(wg, TM, ma);
	for(int i=0; i<=ma.num-1; ++i)
	{
		fprintf(stderr, "[%d] %s: neff = ", i, ma(i).ids); ma(i).neff.nice(stderr); fprintf(stderr, ", Lp = %g\n", ma(i).Lp()); 
		ma(i).plot(HY, REP, disp, 500, dig10(i), dig1(i), 'L');
		ma(i).plot(HY, IMP, disp, 500, dig10(i), dig1(i), 'V');
		ma(i).plot(HY, MOD, disp, 500, dig10(i), dig1(i), 'L');
	}

	fprintf(stderr, "Ok.\n");
	return 0;
}