/* 
 * METRIC-Application:
 * rib waveguide, VEIMS mode analysis 
 */

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

#define Pol   TE     // light polarization
#define Wgpns 1.45   // substrate refractive index
#define Wgpnf 1.65   // guiding film, refractive index
#define Wgpnc 1.0    // cover: air
#define Wgpw  2.5    // rib width
#define Wgph  1.0    // film height
#define Wgpd  0.270  // etching depth
#define Wgpl  1.550  // vacuum wavelength 

/* plots: display window */
Rect Display(-Wgph, -1.1*Wgpw, Wgph*1.5, 1.1*Wgpw);

/* waveguide definition */
WgCrs wgdef()
{
	WgCrs g(2, 1);

	g.hx(0) = 0.0;
	g.hx(1) = Wgph-Wgpd;
	g.hx(2) = Wgph;
	g.hy(0) = -Wgpw/2.0;
	g.hy(1) =  Wgpw/2.0;

	g.n(0,0) = Wgpns;
	g.n(0,1) = Wgpns;
	g.n(0,2) = Wgpns;
	g.n(1,0) = Wgpnf;
	g.n(1,1) = Wgpnf;
	g.n(1,2) = Wgpnf;
	g.n(2,0) = Wgpnc;
	g.n(2,1) = Wgpnf;
	g.n(2,2) = Wgpnc;
	g.n(3,0) = Wgpnc;
	g.n(3,1) = Wgpnc;
	g.n(3,2) = Wgpnc;

	g.lambda = Wgpl;

	return g;
}

/* VEIMS mode analysis, mode profile plots */
int main()
{
	// waveguide definition, refractive index plot
	WgCrs wg = wgdef();
	wg.plot(Display, '0', '0');

	// containers for modes and mode arrays
	EIMode m;	
	EIModeArray ma;

	// invoke the VEIMS mode solver, found modes -> ma 	
	veims(wg, Pol, ma);
	
	for(int j=0; j<=ma.num-1; ++j)
	{
		m = ma(j);
		fprintf(stderr, "[%d]  T%c_(%d, %d): ", 
		        j, polCHR(m.pol), m.vmp.ord, m.hmp.ord);
		fprintf(stderr, "n_eff = %g  ", m.neff);
		fprintf(stderr, "beta = %g\n", m.beta);

		// plots of the constituting 1-D mode profiles
		m.vmp.plot('v', dig1(j)); 
		m.hmp.plot('h', dig1(j)); 

		// mode profile plots, different variants		
		Fcomp cp = principalcomp(m.pol);
		m.plot(cp, ORG, Display, 150, 150, dig10(j), dig1(j), 'C');
		m.plot(cp, SQR, Display, 150, 150, dig10(j), dig1(j), 'C');
		m.plot(cp, ORG, Display, 75, 75, dig10(j), dig1(j), 'S');
		m.plot(cp, ORG, Display, 150, 150, dig10(j), dig1(j), 'I');
		m.fplot(cp, Display, 150, 150, dig10(j), dig1(j));
		m.acplot(Display, 75, 75, dig10(j), dig1(j));

		// generate a Matlab mode profile viewer
		m.viewer(Display, 150, 150, dig10(j), dig1(j));
	}
	fprintf(stderr, "Ok.\n");
	return 0;
}