#include "omm.h"

typedef struct {int xy, i, j, on;} segment_t;

void setgeometry(int nx, int ny, double lx, double ly, double h, double r, int sdiv, int lgeometry, int nseg, segment_t seg[])
{
	assert((nx > 0) && (ny > 0) && (sdiv > 0) && (lgeometry > 0) && (nseg > 0));

	const double x0 = -lx / 2;
	const double y0 = -ly / 2;
	const double dx = lx / nx;
	const double dy = ly / ny;

	// initialize
	for (int n = 0; n < lgeometry; n++) {
		Geometry[n].shape = UNDEF;
		//Geometry[n].div[0] = 1;
		Geometry[n].ifeed = 0;
		Geometry[n].iload = 0;
		//Geometry[n].iradius = 0;
		Geometry[n].iradius = 1;
		Geometry[n].radius = r;
		Geometry[n].offset[0] = Geometry[n].offset[1] = Geometry[n].offset[2] = 0;
	}

	// geometry id
	int gid = 0;

	// feed elements
	const double d = (lx + ly) / (nx + ny);
	const int nefeed = MAX(NINT(sdiv * h, d), 1);

	// feed
	const double xf = 0;
	const double yf = 0;
	Geometry[gid].shape = WIRE;
	Geometry[gid].cosys = XYZ;
	Geometry[gid].pos[0][0] =
	Geometry[gid].pos[1][0] = xf;
	Geometry[gid].pos[0][1] =
	Geometry[gid].pos[1][1] = yf;
	Geometry[gid].pos[0][2] = 0;
	Geometry[gid].pos[1][2] = h / nefeed;
	Geometry[gid].div[0] = 1;
	Geometry[gid].ifeed = 1;
	Geometry[gid].feed[0] = 1;
	Geometry[gid].feed[1] = 0;
	gid++;

	if (nefeed > 1) {
		Geometry[gid].shape = WIRE;
		Geometry[gid].cosys = XYZ;
		Geometry[gid].pos[0][0] =
		Geometry[gid].pos[1][0] = xf;
		Geometry[gid].pos[0][1] =
		Geometry[gid].pos[1][1] = yf;
		Geometry[gid].pos[0][2] = h / nefeed;
		Geometry[gid].pos[1][2] = h;
		Geometry[gid].div[0] = nefeed - 1;
		gid++;
	}

	// patch
	for (int n = 0; n < nseg; n++) {
		if (seg[n].on) {
			const int i = seg[n].i;
			const int j = seg[n].j;
			Geometry[gid].shape = WIRE;
			Geometry[gid].cosys = XYZ;
			assert((seg[n].xy == 1) || (seg[n].xy == 2));
			if      (seg[n].xy == 1) {
				assert(i >= 0 && i < nx && j >= 0 && j <= ny);
				Geometry[gid].pos[0][0] = x0 + (i + 0) * dx;
				Geometry[gid].pos[1][0] = x0 + (i + 1) * dx;
				Geometry[gid].pos[0][1] =
				Geometry[gid].pos[1][1] = y0 + j * dy;
			}
			else if (seg[n].xy == 2) {
				assert(i >= 0 && i <= nx && j >= 0 && j < ny);
				Geometry[gid].pos[0][0] =
				Geometry[gid].pos[1][0] = x0 + i * dx;
				Geometry[gid].pos[0][1] = y0 + (j + 0) * dy;
				Geometry[gid].pos[1][1] = y0 + (j + 1) * dy;
			}
			Geometry[gid].pos[0][2] =
			Geometry[gid].pos[1][2] = h;
			Geometry[gid].div[0] = sdiv;
			gid++;
		}
	}

	NGeometry = gid;
}
