/*
TDS.c
テーパースロットアンテナ
*/

#include <stdio.h>
#include <math.h>
#include "ofd_datalib.h"

static void fa(double, double, double, double, double, double, double, double, double, double *, double *);
static void fb(double, double, double, double, double, double, double, double, double, double *, double *);
static void fc(double, double, double, double, double, double, double, double, double, double, double *, double *);

int main(void)
{
	const char AntType = 'A';    // アンテナタイプ: 'A'/'B'/'C'
	const double W = 200e-3;
	const double L = 220e-3;
	const double TW = 40e-3;
	const double SW = 1.2e-3;
	const double TL = 24e-3;
	const double TH = 26e-3;
	const double FL = 56e-3;
	const double margin = 50e-3;  // 余白
	const double dx = 1.2e-3;     // X方向基本セルサイズ
	const double dy = 1.2e-3;     // Y方向基本セルサイズ
	const double dz = 1.2e-3;     // Z方向基本セルサイズ

	char title[BUFSIZ], outfile[BUFSIZ];
	sprintf(title, "taper slot antenna: Type %c", AntType);
	sprintf(outfile, "TSA_%c.ofd", AntType);

	// initialize

	ofd_init();

	// title

	ofd_title(title);

	// mesh

	const double x0 = - W / 2 - margin;
	const double x1 = - W / 2;
	const double x2 = - SW / 2;
	const double x3 = + SW / 2;
	const double x4 = + W / 2;
	const double x5 = + W / 2 + margin;
	ofd_xsection(6, x0, x1, x2, x3, x4, x5);
	ofd_xdivision(5, NINT(x1 - x0, dx), NINT(x2 - x1, dx), 1, NINT(x4 - x3, dx), NINT(x5 - x4, dx));

	const double y0 = - margin;
	const double y1 = 0;
	const double y2 = L;
	const double y3 = L + margin;
	ofd_ysection(4, y0, y1, y2, y3);
	ofd_ydivision(3, NINT(y1 - y0, dy), NINT(y2 - y1, dy), NINT(y3 - y2, dy));

	const double z0 = - margin;
	const double z1 = 0;
	const double z2 = + margin;
	ofd_zsection(3, z0, z1, z2);
	ofd_zdivision(2, NINT(z1 - z0, dz), NINT(z2 - z1, dz));

	// material

	//ofd_material(2.0, 0.0, 1.0, 0.0, "");

	// geometry

	double L1 = 0, L2 = 0, W1 = 0, W2 = 0;
	if      (AntType == 'A') {
		L1 = 140e-3;
		W1 = (W - SW) / 2;
	}
	else if (AntType == 'B') {
		L1 = 140e-3;
		W1 = (W - SW) / 2;
	}
	else if (AntType == 'C') {
		L1 = 40e-3;
		L2 = (L - 4 * L1) / 2;
		W2 = 35e-3;
		W1 = (W - SW) / 2 - W2;
	}

	const int m = 1;
	const int shape = 1;
	const double gz1 = 0e-3;
	const double gz2 = 0e-3;
	const int ny = NINT(L, dy);
	for (int j = 0; j < ny; j++) {
		double X1 = 0, X2 = 0;
		const double gy1 = (j + 0) * dy;
		const double gy2 = (j + 1) * dy;
		const double y = (gy1 + gy2) / 2;
		if      (AntType == 'A') {
			fa(y, W, TW, SW, TL, TH, L, L1, W1, &X1, &X2);
		}
		else if (AntType == 'B') {
			fb(y, W, TW, SW, TL, TH, L, L1, W1, &X1, &X2);
		}
		else if (AntType == 'C') {
			fc(y, W, TW, SW, TL, TH, L1, W1, L2, W2, &X1, &X2);
		}
		ofd_geometry(m, shape, +X1, +X2, gy1, gy2, gz1, gz2);
		ofd_geometry(m, shape, -X1, -X2, gy1, gy2, gz1, gz2);
	}

	// feed

	ofd_feed('X', 0e-3, FL, 0e-3, 1, 0, 50);
	//ofd_rfeed(10);

	// ABC

	//ofd_pml(5, 2, 1e-5);

	// frequency

	ofd_frequency1(0.5e9, 4.5e9, 20);
	ofd_frequency2(1e9, 2e9, 1);

	// solver

	ofd_solver(10000, 100, 1e-3);

	// iteration

	ofd_plotiter(1);

	// waveform and spectrum

	//ofd_plotfeed(1);
	//ofd_plotpoint(1);

	// frequency

	//ofd_plotsmith(1);
	//ofd_plotzin(1, 0, 0, 0);
	//ofd_plotyin(1, 0, 0, 0);
	ofd_plotref(2, -30, 0, 6);
	//ofd_plotfar0d(90, 0, 1, 0, 0, 0);
	ofd_freqdiv(8);

	// far-1d

	ofd_plotfar1d('Z', 180, 0);
	ofd_far1dstyle(0);
	ofd_far1dcomponent(1, 0, 0);
	ofd_far1ddb(1);
	ofd_far1dnorm(0);
	ofd_far1dabs(1);
	ofd_far1dscale(-20, +10, 6);

	// far-2d

	//ofd_plotfar2d(18, 36);
	//ofd_far2dcomponent(1, 0, 0, 0, 0, 0, 0);
	//ofd_far2ddb(1);
	//ofd_far2dscale(-30, +10);
	//ofd_far2dobj(0.5);

	// near-1d

	//ofd_plotnear1d("E", 'Z', 30e-3, 0e-3);
	//ofd_near1ddb(1);
	//ofd_near1dscale(-30, +10, 4);

	// near-2d

	//ofd_plotnear2d("E", 'X', 30e-3);
	//ofd_near2ddim(1, 1);
	//ofd_near2ddb(1);
	//ofd_near2dscale(-30, +10);
	//ofd_near2dcontour(0);
	//ofd_near2dobj(1);
	//ofd_near2dzoom(-50e-3, 50e-3, -50e-3, 50e-3);

	// output

	ofd_outdata(outfile);

	return 0;
}


// アンテナ形状関数(Type A)
static void fa(double y, double W, double TW, double SW, double TL, double TH, double L, double L1, double W1, double *X1, double *X2)
{
	if      (y < TL) {
		*X1 = 0e-3;
	}
	else if (y < TL + TH) {
		*X1 = TW / 2 - (y - TL) / TH * (TW - SW) / 2;
	}
	else if (y < L - L1) {
		*X1 = SW / 2;
	}
	else {
		*X1 = W / 2 - (L - y) / L1 * W1;
	}

	*X2 = W / 2;
}


// アンテナ形状関数(Type B)
static void fb(double y, double W, double TW, double SW, double TL, double TH, double L, double L1, double W1, double *X1, double *X2)
{
	if      (y < TL) {
		*X1 = 0e-3;
	}
	else if (y < TL + TH) {
		*X1 = TW / 2 - (y - TL) / TH * (TW - SW) / 2;
	}
	else if (y < L - L1) {
		*X1 = SW / 2;
	}
	else {
		const double fy = (y - L + L1) / L1;
		*X1 = W / 2 - W1 * sqrt(1 - fy * fy);
	}

	*X2 = W / 2;
}


//アンテナ形状関数(Type C)
static void fc(double y, double W, double TW, double SW, double TL, double TH, double L1, double W1, double L2, double W2, double *X1, double *X2)
{
	if      (y < TL) {
		*X1 = 0e-3;
	}
	else if (y < TL + TH) {
		*X1 = TW / 2 - (y - TL) / TH * (TW - SW) / 2;
	}
	else if (y < 2 * L1) {
		*X1 = SW / 2;
	}
	else {
		const double fy = (y - 2 * L1) / (2 * (L1 + L2));
		*X1 = (W / 2 - W2) - W1 * sqrt(1 - fy * fy);
	}

	if      (y < 2 * L1) {
		const double fy = (y - L1) / L1;
		*X2 = (W / 2 - W2) + W2 * sqrt(1 - fy * fy);
	}
	else if (y < 2 * (L1 + L2)) {
		const double fy = (y - 2 * L1 - L2) / L2;
		*X2 = (W / 2 - W2) - W2 * sqrt(1 - fy * fy);
	}
	else {
		const double fy = (y - 3 * L1 - 2 * L2) / L1;
		*X2 = (W / 2 - W2) + W2 * sqrt(1 - fy * fy);
	}
}