diff --git a/src/context/minimal.rs b/src/context/minimal.rs
index 194be4f..9b6d485 100644
--- a/src/context/minimal.rs
+++ b/src/context/minimal.rs
@@ -136,11 +136,17 @@ mod tests {
         ctx.register_resource("stupid:way", "addone | addone | addone inv");
         let op = ctx.op("stupid:way")?;
 
+        let steps = ctx.steps(op)?;
+        assert_eq!(steps.len(), 3);
+        assert_eq!(steps[0], "addone");
+        assert_eq!(steps[1], "addone");
+        assert_eq!(steps[2], "addone inv");
+
         let mut data = some_basic_coor2dinates();
         assert_eq!(data[0][0], 55.);
         assert_eq!(data[1][0], 59.);
 
-        ctx.apply(op, Fwd, &mut data)?;
+        assert_eq!(2, ctx.apply(op, Fwd, &mut data)?);
         assert_eq!(data[0][0], 56.);
         assert_eq!(data[1][0], 60.);
 
@@ -148,12 +154,6 @@ mod tests {
         assert_eq!(data[0][0], 55.);
         assert_eq!(data[1][0], 59.);
 
-        let steps = ctx.steps(op)?;
-        assert_eq!(steps.len(), 3);
-        assert_eq!(steps[0], "addone");
-        assert_eq!(steps[1], "addone");
-        assert_eq!(steps[2], "addone inv");
-
         let params = ctx.params(op, 1)?;
         let ellps = params.ellps(0);
         assert_eq!(ellps.semimajor_axis(), 6378137.);
@@ -195,9 +195,9 @@ mod tests {
         // while ellps_? defaults to GRS80 - so they are there even though we havent
         // set them
         let params = ctx.params(op, 1)?;
-        let ellps = params.ellps[0];
+        let ellps = params.ellps(0);
         assert_eq!(ellps.semimajor_axis(), 6378137.);
-        assert_eq!(0., params.lat[0]);
+        assert_eq!(0., params.real("lat_0")?);
 
         // The zone id is found among the natural numbers (which here includes 0)
         let zone = params.natural("zone")?;
@@ -225,7 +225,7 @@ mod tests {
         let op = ctx.op("utm zone=32")?;
         let steps = ctx.steps(op)?;
         assert!(steps.len() == 1);
-        let ellps = ctx.params(op, 0)?.ellps[0];
+        let ellps = ctx.params(op, 0)?.ellps(0);
         let jac = Jacobian::new(
             &ctx,
             op,
diff --git a/src/inner_op/btmerc.rs b/src/inner_op/btmerc.rs
index 0d4ca47..98882b0 100644
--- a/src/inner_op/btmerc.rs
+++ b/src/inner_op/btmerc.rs
@@ -5,13 +5,13 @@ use crate::operator_authoring::*;
 
 // Forward transverse mercator, following Bowring (1989)
 fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
-    let ellps = op.params.ellps[0];
+    let ellps = op.params.ellps(0);
     let eps = ellps.second_eccentricity_squared();
-    let lat_0 = op.params.lat[0];
-    let lon_0 = op.params.lon[0];
-    let x_0 = op.params.x[0];
-    let y_0 = op.params.y[0];
-    let k_0 = op.params.k[0];
+    let lat_0 = op.params.lat(0).to_radians();
+    let lon_0 = op.params.lon(0).to_radians();
+    let x_0 = op.params.x(0);
+    let y_0 = op.params.y(0);
+    let k_0 = op.params.k(0);
 
     let mut successes = 0_usize;
     let n = operands.len();
@@ -53,13 +53,13 @@ fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
 
 // Inverse transverse mercator, following Bowring (1989)
 fn inv(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
-    let ellps = op.params.ellps[0];
+    let ellps = op.params.ellps(0);
     let eps = ellps.second_eccentricity_squared();
-    let lat_0 = op.params.lat[0];
-    let lon_0 = op.params.lon[0];
-    let x_0 = op.params.x[0];
-    let y_0 = op.params.y[0];
-    let k_0 = op.params.k[0];
+    let lat_0 = op.params.lat(0).to_radians();
+    let lon_0 = op.params.lon(0).to_radians();
+    let x_0 = op.params.x(0);
+    let y_0 = op.params.y(0);
+    let k_0 = op.params.k(0);
 
     let mut successes = 0_usize;
     let n = operands.len();
@@ -134,22 +134,24 @@ pub fn utm(parameters: &RawParameters, _ctx: &dyn Context) -> Result<Op, Error>
     }
 
     // The scaling factor is 0.9996 by definition of UTM
-    params.k[0] = 0.9996;
+    params.real.insert("k_0", 0.9996);
 
     // The center meridian is determined by the zone
-    params.lon[0] = (-183. + 6. * zone as f64).to_radians();
+    params
+        .real
+        .insert("lon_0", -183. + 6. * zone as f64);
 
     // The base parallel is by definition the equator
-    params.lat[0] = 0.0;
+    params.real.insert("lat_0", 0.);
 
     // The false easting is 500000 m by definition of UTM
-    params.x[0] = 500000.0;
+    params.real.insert("x_0", 500_000.);
 
     // The false northing is 0 m by definition of UTM
-    params.y[0] = 0.0;
+    params.real.insert("y_0", 0.);
     // or 10_000_000 m if using the southern aspect
     if params.boolean("south") {
-        params.y[0] = 10_000_000.0;
+        params.real.insert("y_0", 10_000_000.0);
     }
 
     let descriptor = OpDescriptor::new(def, InnerOp(fwd), Some(InnerOp(inv)));
diff --git a/src/inner_op/cart.rs b/src/inner_op/cart.rs
index b9fe5ed..bf421fe 100644
--- a/src/inner_op/cart.rs
+++ b/src/inner_op/cart.rs
@@ -6,9 +6,10 @@ use crate::operator_authoring::*;
 fn cart_fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
     let n = operands.len();
     let mut successes = 0;
+    let ellps = op.params.ellps(0);
     for i in 0..n {
         let mut coord = operands.get_coord(i);
-        coord = op.params.ellps[0].cartesian(&coord);
+        coord = ellps.cartesian(&coord);
         if !coord.0.iter().any(|c| c.is_nan()) {
             successes += 1;
         }
@@ -20,20 +21,21 @@ fn cart_fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> us
 // ----- I N V E R S E --------------------------------------------------------------
 
 fn cart_inv(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
+    let ellps = op.params.ellps(0);
     // eccentricity squared, Fukushima's E, Claessens' c3 = 1-c2`
-    let es = op.params.ellps[0].eccentricity_squared();
+    let es = ellps.eccentricity_squared();
     // semiminor axis
-    let b = op.params.ellps[0].semiminor_axis();
+    let b = ellps.semiminor_axis();
     // semimajor axis
-    let a = op.params.ellps[0].semimajor_axis();
+    let a = ellps.semimajor_axis();
     // reciproque of a
-    let ra = 1. / op.params.ellps[0].semimajor_axis();
+    let ra = 1. / ellps.semimajor_axis();
     // aspect ratio, b/a: Fukushima's ec, Claessens' c4
     let ar = b * ra;
     // 1.5 times the fourth power of the eccentricity
     let ce4 = 1.5 * es * es;
     // if we're closer than this to the Z axis, we force latitude to one of the poles
-    let cutoff = op.params.ellps[0].semimajor_axis() * 1e-16;
+    let cutoff = ellps.semimajor_axis() * 1e-16;
 
     let n = operands.len();
     let mut successes = 0;
diff --git a/src/inner_op/curvature.rs b/src/inner_op/curvature.rs
index a87443f..691165e 100644
--- a/src/inner_op/curvature.rs
+++ b/src/inner_op/curvature.rs
@@ -6,7 +6,7 @@ use crate::operator_authoring::*;
 fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
     let n = operands.len();
     let sliced = 0..n;
-    let ellps = op.params.ellps[0];
+    let ellps = op.params.ellps(0);
 
     let prime = op.params.boolean("prime");
     let meridional = op.params.boolean("meridian");
@@ -116,6 +116,11 @@ pub fn new(parameters: &RawParameters, ctx: &dyn Context) -> Result<Op, Error> {
         ));
     }
 
+    // Check that we have a proper ellipsoid
+    if op.params.text("ellps").is_ok() {
+        let _ = Ellipsoid::named(op.params.text("ellps")?.as_str())?;
+    }
+
     Ok(op)
 }
 
diff --git a/src/inner_op/geodesic.rs b/src/inner_op/geodesic.rs
index d218b60..3e0e9d6 100644
--- a/src/inner_op/geodesic.rs
+++ b/src/inner_op/geodesic.rs
@@ -4,7 +4,7 @@ use crate::operator_authoring::*;
 // ----- F O R W A R D -----------------------------------------------------------------
 
 fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
-    let ellps = op.params.ellps[0];
+    let ellps = op.params.ellps(0);
 
     let n = operands.len();
     let sliced = 0..n;
@@ -35,7 +35,7 @@ fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
 // ----- I N V E R S E -----------------------------------------------------------------
 
 fn inv(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
-    let ellps = op.params.ellps[0];
+    let ellps = op.params.ellps(0);
     let reversible = op.params.boolean("reversible");
 
     let n = operands.len();
diff --git a/src/inner_op/laea.rs b/src/inner_op/laea.rs
index c530896..47f7cf9 100644
--- a/src/inner_op/laea.rs
+++ b/src/inner_op/laea.rs
@@ -20,9 +20,9 @@ fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
     let north_polar = op.params.boolean("north_polar");
     let south_polar = op.params.boolean("south_polar");
 
-    let lon_0 = op.params.lon(0);
-    let x_0 = op.params.x(0);
-    let y_0 = op.params.y(0);
+    let lon_0 = op.params.real("lon_0").unwrap_or(0.).to_radians();
+    let x_0 = op.params.real("x_0").unwrap_or(0.);
+    let y_0 = op.params.real("y_0").unwrap_or(0.);
     let ellps = op.params.ellps(0);
     let e = ellps.eccentricity();
     let a = ellps.semimajor_axis();
@@ -95,10 +95,10 @@ fn inv(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
     let north_polar = op.params.boolean("north_polar");
     let south_polar = op.params.boolean("south_polar");
 
-    let lon_0 = op.params.lon(0);
-    let lat_0 = op.params.lat(0);
-    let x_0 = op.params.x(0);
-    let y_0 = op.params.y(0);
+    let lon_0 = op.params.real("lon_0").unwrap_or(0.).to_radians();
+    let lat_0 = op.params.real("lat_0").unwrap_or(0.).to_radians();
+    let x_0 = op.params.real("x_0").unwrap_or(0.);
+    let y_0 = op.params.real("y_0").unwrap_or(0.);
 
     let ellps = op.params.ellps(0);
     let a = ellps.semimajor_axis();
@@ -191,7 +191,8 @@ pub fn new(parameters: &RawParameters, _ctx: &dyn Context) -> Result<Op, Error>
     let def = &parameters.definition;
     let mut params = ParsedParameters::new(parameters, &GAMUT)?;
 
-    let lat_0 = params.lat[0];
+    let lat_0 = params.real("lat_0").unwrap_or(0.).to_radians();
+
     if lat_0.is_nan() {
         warn!("LAEA: Bad central latitude!");
         return Err(Error::BadParam("lat_0".to_string(), def.clone()));
@@ -216,7 +217,7 @@ pub fn new(parameters: &RawParameters, _ctx: &dyn Context) -> Result<Op, Error>
 
     // --- Precompute some latitude invariant factors ---
 
-    let ellps = params.ellps[0];
+    let ellps = params.ellps(0);
     let a = ellps.semimajor_axis();
     let es = ellps.eccentricity_squared();
     let e = es.sqrt();
@@ -266,6 +267,7 @@ pub fn new(parameters: &RawParameters, _ctx: &dyn Context) -> Result<Op, Error>
 #[cfg(test)]
 mod tests {
     use super::*;
+    use float_eq::assert_float_eq;
 
     #[test]
     fn laea_oblique() -> Result<(), Error> {
@@ -275,11 +277,16 @@ mod tests {
         let op = ctx.op("laea ellps=GRS80 lat_0=52 lon_0=10  x_0=4321000 y_0=3210000")?;
 
         // The test point from IOGP
-        let p = Coor4D::geo(50.0, 5.0, 0.0, 0.0);
-        let mut operands = [p];
+        let p = Coor2D::geo(50.0, 5.0);
+        let geo = [p];
+        let p = Coor2D::raw(3962799.45, 2999718.85);
+        let projected = [p];
+
+        let mut operands = geo.clone();
 
         // Forward
         ctx.apply(op, Fwd, &mut operands)?;
+        assert_float_eq!(operands[0].0, projected[0].0, abs_all <= 0.01);
         assert!((operands[0][0] - 3962799.45).abs() < 0.01);
         assert!((operands[0][1] - 2999718.85).abs() < 0.01);
         ctx.apply(op, Inv, &mut operands)?;
diff --git a/src/inner_op/latitude.rs b/src/inner_op/latitude.rs
index 8f324e3..16f007c 100644
--- a/src/inner_op/latitude.rs
+++ b/src/inner_op/latitude.rs
@@ -7,7 +7,7 @@ fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
     let mut successes = 0_usize;
     let n = operands.len();
     let sliced = 0..n;
-    let ellps = op.params.ellps[0];
+    let ellps = op.params.ellps(0);
 
     if op.params.boolean("geocentric") {
         for i in sliced {
@@ -66,7 +66,7 @@ fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
 fn inv(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
     let mut successes = 0_usize;
     let n = operands.len();
-    let ellps = op.params.ellps[0];
+    let ellps = op.params.ellps(0);
 
     if op.params.boolean("geocentric") {
         for i in 0..n {
@@ -137,7 +137,7 @@ pub const GAMUT: [OpParameter; 8] = [
 
 pub fn new(parameters: &RawParameters, ctx: &dyn Context) -> Result<Op, Error> {
     let mut op = Op::plain(parameters, InnerOp(fwd), Some(InnerOp(inv)), &GAMUT, ctx)?;
-    let ellps = op.params.ellps[0];
+    let ellps = op.params.ellps(0);
 
     let mut number_of_flags = 0_usize;
     if op.params.boolean("geocentric") {
diff --git a/src/inner_op/lcc.rs b/src/inner_op/lcc.rs
index 5e4e49b..1eed9df 100644
--- a/src/inner_op/lcc.rs
+++ b/src/inner_op/lcc.rs
@@ -9,12 +9,13 @@ const EPS10: f64 = 1e-10;
 // Forward Lambert conformal conic, following the PROJ implementation,
 // cf.  https://proj.org/operations/projections/lcc.html
 fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
-    let a = op.params.ellps[0].semimajor_axis();
-    let e = op.params.ellps[0].eccentricity();
-    let lon_0 = op.params.lon[0];
-    let k_0 = op.params.k[0];
-    let x_0 = op.params.x[0];
-    let y_0 = op.params.y[0];
+    let ellps = op.params.ellps(0);
+    let a = ellps.semimajor_axis();
+    let e = ellps.eccentricity();
+    let lon_0 = op.params.lon(0);
+    let k_0 = op.params.k(0);
+    let x_0 = op.params.x(0);
+    let y_0 = op.params.y(0);
     let Ok(n) = op.params.real("n") else { return 0 };
     let Ok(c) = op.params.real("c") else { return 0 };
     let Ok(rho0) = op.params.real("rho0") else { return 0 };
@@ -48,12 +49,13 @@ fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
 
 // ----- I N V E R S E -----------------------------------------------------------------
 fn inv(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
-    let a = op.params.ellps[0].semimajor_axis();
-    let e = op.params.ellps[0].eccentricity();
-    let lon_0 = op.params.lon[0];
-    let k_0 = op.params.k[0];
-    let x_0 = op.params.x[0];
-    let y_0 = op.params.y[0];
+    let ellps = op.params.ellps(0);
+    let a = ellps.semimajor_axis();
+    let e = ellps.eccentricity();
+    let lon_0 = op.params.lon(0);
+    let k_0 = op.params.k(0);
+    let x_0 = op.params.x(0);
+    let y_0 = op.params.y(0);
     let Ok(n) = op.params.real("n") else { return 0 };
     let Ok(c) = op.params.real("c") else { return 0 };
     let Ok(rho0) = op.params.real("rho0") else { return 0 };
@@ -119,18 +121,25 @@ pub const GAMUT: [OpParameter; 9] = [
 pub fn new(parameters: &RawParameters, _ctx: &dyn Context) -> Result<Op, Error> {
     let def = &parameters.definition;
     let mut params = ParsedParameters::new(parameters, &GAMUT)?;
-    if params.lat[2].is_nan() {
-        params.lat[2] = params.lat[1];
+    if !params.real.contains_key("lat_2") {
+        params.real.insert("lat_2", params.lat(1));
     }
 
-    let phi1 = params.lat[1];
-    let mut phi2 = params.lat[2];
+    let phi1 = params.lat(1).to_radians();
+    let mut phi2 = params.lat(2).to_radians();
     if phi2.is_nan() {
         phi2 = phi1;
     }
-    params.lat[2] = phi2;
-
-    let mut lat_0 = params.lat[0];
+    params
+        .real
+        .insert("lon_0", params.real["lon_0"].to_radians());
+    params
+        .real
+        .insert("lat_0", params.real["lat_0"].to_radians());
+    params.real.insert("lat_1", phi1);
+    params.real.insert("lat_2", phi2);
+
+    let mut lat_0 = params.lat(0);
     if lat_0.is_nan() {
         lat_0 = 0.;
         if (phi1 - phi2).abs() < EPS10 {
@@ -140,8 +149,9 @@ pub fn new(parameters: &RawParameters, _ctx: &dyn Context) -> Result<Op, Error>
 
     let sc = phi1.sin_cos();
     let mut n = sc.0;
-    let e = params.ellps[0].eccentricity();
-    let es = params.ellps[0].eccentricity_squared();
+    let ellps = params.ellps(0);
+    let e = ellps.eccentricity();
+    let es = ellps.eccentricity_squared();
 
     if (phi1 + phi2).abs() < EPS10 {
         return Err(Error::General(
@@ -189,7 +199,7 @@ pub fn new(parameters: &RawParameters, _ctx: &dyn Context) -> Result<Op, Error>
     params.real.insert("c", c);
     params.real.insert("n", n);
     params.real.insert("rho0", rho0);
-    params.lat[0] = lat_0;
+    params.real.insert("lat_0", lat_0);
 
     let descriptor = OpDescriptor::new(def, InnerOp(fwd), Some(InnerOp(inv)));
     let steps = Vec::<Op>::new();
diff --git a/src/inner_op/merc.rs b/src/inner_op/merc.rs
index d069ede..69e1372 100644
--- a/src/inner_op/merc.rs
+++ b/src/inner_op/merc.rs
@@ -4,13 +4,13 @@ use crate::operator_authoring::*;
 // ----- F O R W A R D -----------------------------------------------------------------
 
 fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
-    let ellps = op.params.ellps[0];
+    let ellps = op.params.ellps(0);
     let a = ellps.semimajor_axis();
-    let k_0 = op.params.k[0];
-    let x_0 = op.params.x[0];
-    let y_0 = op.params.y[0];
-    let lat_0 = op.params.lat[0];
-    let lon_0 = op.params.lon[0];
+    let k_0 = op.params.k(0);
+    let x_0 = op.params.x(0);
+    let y_0 = op.params.y(0);
+    let lat_0 = op.params.lat(0);
+    let lon_0 = op.params.lon(0);
 
     let mut successes = 0_usize;
     let length = operands.len();
@@ -33,13 +33,13 @@ fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
 // ----- I N V E R S E -----------------------------------------------------------------
 
 fn inv(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
-    let ellps = op.params.ellps[0];
+    let ellps = op.params.ellps(0);
     let a = ellps.semimajor_axis();
-    let k_0 = op.params.k[0];
-    let x_0 = op.params.x[0];
-    let y_0 = op.params.y[0];
-    let lat_0 = op.params.lat[0];
-    let lon_0 = op.params.lon[0];
+    let k_0 = op.params.k(0);
+    let x_0 = op.params.x(0);
+    let y_0 = op.params.y(0);
+    let lat_0 = op.params.lat(0);
+    let lon_0 = op.params.lon(0);
 
     let mut successes = 0_usize;
     let length = operands.len();
@@ -80,7 +80,7 @@ pub const GAMUT: [OpParameter; 8] = [
 pub fn new(parameters: &RawParameters, _ctx: &dyn Context) -> Result<Op, Error> {
     let def = &parameters.definition;
     let mut params = ParsedParameters::new(parameters, &GAMUT)?;
-    let ellps = params.ellps[0];
+    let ellps = params.ellps(0);
 
     let lat_ts = params.real("lat_ts")?;
     if lat_ts.abs() > 90. {
@@ -92,7 +92,8 @@ pub fn new(parameters: &RawParameters, _ctx: &dyn Context) -> Result<Op, Error>
     // lat_ts trumps k_0
     if lat_ts != 0.0 {
         let sc = lat_ts.to_radians().sin_cos();
-        params.k[0] = sc.1 / (1. - ellps.eccentricity_squared() * sc.0 * sc.0).sqrt()
+        let k_0 = sc.1 / (1. - ellps.eccentricity_squared() * sc.0 * sc.0).sqrt();
+        params.real.insert("k_0", k_0);
     }
 
     let descriptor = OpDescriptor::new(def, InnerOp(fwd), Some(InnerOp(inv)));
diff --git a/src/inner_op/molodensky.rs b/src/inner_op/molodensky.rs
index 35af18f..b6c8272 100644
--- a/src/inner_op/molodensky.rs
+++ b/src/inner_op/molodensky.rs
@@ -25,7 +25,7 @@ fn common(
     operands: &mut dyn CoordinateSet,
     direction: Direction,
 ) -> usize {
-    let ellps = op.params.ellps[0];
+    let ellps = op.params.ellps(0);
     let a = ellps.semimajor_axis();
     let f = ellps.flattening();
     let es = ellps.eccentricity_squared();
@@ -100,8 +100,8 @@ pub fn new(parameters: &RawParameters, _ctx: &dyn Context) -> Result<Op, Error>
     let def = &parameters.definition;
     let mut params = ParsedParameters::new(parameters, &GAMUT)?;
 
-    let ellps_0 = params.ellps[0];
-    let ellps_1 = params.ellps[1];
+    let ellps_0 = params.ellps(0);
+    let ellps_1 = params.ellps(1);
 
     // We may use `ellps, da, df`, to parameterize the op, but `ellps_0, ellps_1`
     // is a more likely set of parameters to come across in real life.
diff --git a/src/inner_op/omerc.rs b/src/inner_op/omerc.rs
index 0355dee..1929ce8 100644
--- a/src/inner_op/omerc.rs
+++ b/src/inner_op/omerc.rs
@@ -9,14 +9,14 @@ use std::f64::consts::FRAC_PI_4;
 
 #[allow(non_snake_case)]
 fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
-    let ellps = op.params.ellps[0];
+    let ellps = op.params.ellps(0);
     let es = ellps.eccentricity_squared();
     let e = es.sqrt();
 
-    let kc = op.params.k[0];
+    let kc = op.params.k(0);
 
-    let FE = op.params.x[0];
-    let FN = op.params.y[0];
+    let FE = op.params.x(0);
+    let FN = op.params.y(0);
     let Ec = FE;
     let Nc = FN;
 
@@ -92,6 +92,7 @@ fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
             let u = A * (S * c0 + V * s0).atan2(cblon) / B;
             coord[0] = v * cc + u * sc + FE;
             coord[1] = u * cc - v * sc + FN;
+            successes += 1;
             continue;
         }
 
@@ -106,6 +107,7 @@ fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
             };
             coord[0] = v * cc + u * sc + Ec;
             coord[1] = u * cc - v * sc + Nc;
+            successes += 1;
             continue;
         }
 
@@ -125,14 +127,14 @@ fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
 
 #[allow(non_snake_case)]
 fn inv(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
-    let ellps = op.params.ellps[0];
+    let ellps = op.params.ellps(0);
     let es = ellps.eccentricity_squared();
     let e = es.sqrt();
 
-    let kc = op.params.k[0];
+    let kc = op.params.k(0);
 
-    let FE = op.params.x[0];
-    let FN = op.params.y[0];
+    let FE = op.params.x(0);
+    let FN = op.params.y(0);
 
     let latc = op.params.real["latc"].to_radians();
     let lonc = op.params.real["lonc"].to_radians();
@@ -276,6 +278,7 @@ pub fn new(parameters: &RawParameters, _ctx: &dyn Context) -> Result<Op, Error>
 #[cfg(test)]
 mod tests {
     use super::*;
+    use float_eq::assert_float_eq;
 
     #[test]
     fn omerc() -> Result<(), Error> {
@@ -291,15 +294,27 @@ mod tests {
 
         // Validation value from EPSG
         let geo = [Coor2D::geo(5.3872535833, 115.8055054444)];
-
         let projected = [Coor2D::raw(679245.7281740266, 596562.7774687681)];
 
         // Forward
         let mut operands = geo.clone();
 
+        assert_eq!(1, ctx.apply(op, Fwd, &mut operands)?);
+        for i in 0..operands.len() {
+            assert_float_eq!(operands[i].0, projected[i].0, abs_all <= 1e-9);
+        }
+
+        // Roundtrip
+        assert_eq!(1, ctx.apply(op, Inv, &mut operands)?);
+        for i in 0..operands.len() {
+            assert_float_eq!(operands[i].0, geo[i].0, abs_all <= 1e-9);
+        }
+
+        // Forward
+        let mut operands = geo.clone();
+
         ctx.apply(op, Fwd, &mut operands)?;
         for i in 0..operands.len() {
-            dbg!(operands[i]);
             assert!(operands[i].hypot2(&projected[i]) < 1e-9);
         }
 
diff --git a/src/inner_op/pipeline.rs b/src/inner_op/pipeline.rs
index 8f3746e..2b8558b 100644
--- a/src/inner_op/pipeline.rs
+++ b/src/inner_op/pipeline.rs
@@ -194,7 +194,7 @@ mod tests {
         let op = ctx.op("addone|addone|addone")?;
         let mut data = some_basic_coor2dinates();
 
-        ctx.apply(op, Fwd, &mut data)?;
+        assert_eq!(2, ctx.apply(op, Fwd, &mut data)?);
         assert_eq!(data[0][0], 58.);
         assert_eq!(data[1][0], 62.);
 
diff --git a/src/inner_op/tmerc.rs b/src/inner_op/tmerc.rs
index f660f4a..71543b9 100644
--- a/src/inner_op/tmerc.rs
+++ b/src/inner_op/tmerc.rs
@@ -7,10 +7,10 @@ use crate::operator_authoring::*;
 // Forward transverse mercator, following Engsager & Poder(2007)
 fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
     // Make all precomputed parameters directly accessible
-    let ellps = op.params.ellps[0];
-    let lat_0 = op.params.lat[0];
-    let lon_0 = op.params.lon[0];
-    let x_0 = op.params.x[0];
+    let ellps = op.params.ellps(0);
+    let lat_0 = op.params.lat(0).to_radians();
+    let lon_0 = op.params.lon(0).to_radians();
+    let x_0 = op.params.x(0);
     let Some(conformal) = op.params.fourier_coefficients.get("conformal") else {
         warn!("Missing Fourier coefficients for conformal mapping!");
         return 0;
@@ -94,12 +94,12 @@ fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
 
 // ----- I N V E R S E -----------------------------------------------------------------
 
-// Inverse Transverse Mercator, following Engsager & Poder (2007) (currently Bowring stands in!)
+// Inverse Transverse Mercator, following Engsager & Poder (2007)
 fn inv(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
     // Make all precomputed parameters directly accessible
-    let ellps = op.params.ellps[0];
-    let lon_0 = op.params.lon[0];
-    let x_0 = op.params.x[0];
+    let ellps = op.params.ellps(0);
+    let lon_0 = op.params.lon(0).to_radians();
+    let x_0 = op.params.x(0);
     let Some(conformal) = op.params.fourier_coefficients.get("conformal") else {
         warn!("Missing Fourier coefficients for conformal mapping!");
         return 0;
@@ -201,22 +201,24 @@ pub fn utm(parameters: &RawParameters, _ctx: &dyn Context) -> Result<Op, Error>
     }
 
     // The scaling factor is 0.9996 by definition of UTM
-    params.k[0] = 0.9996;
+    params.real.insert("k_0", 0.9996);
 
     // The center meridian is determined by the zone
-    params.lon[0] = (-183. + 6. * zone as f64).to_radians();
+    params
+        .real
+        .insert("lon_0", -183. + 6. * zone as f64);
 
     // The base parallel is by definition the equator
-    params.lat[0] = 0.0;
+    params.real.insert("lat_0", 0.);
 
     // The false easting is 500000 m by definition of UTM
-    params.x[0] = 500_000.0;
+    params.real.insert("x_0", 500_000.);
 
     // The false northing is 0 m by definition of UTM
-    params.y[0] = 0.0;
+    params.real.insert("y_0", 0.);
     // or 10_000_000 m if using the southern aspect
     if params.boolean("south") {
-        params.y[0] = 10_000_000.0;
+        params.real.insert("y_0", 10_000_000.0);
     }
 
     let descriptor = OpDescriptor::new(def, InnerOp(fwd), Some(InnerOp(inv)));
@@ -263,13 +265,13 @@ const TRANSVERSE_MERCATOR: PolynomialCoefficients = PolynomialCoefficients {
 // Pre-compute some of the computationally heavy prerequisites,
 // to get better amortization over the full operator lifetime.
 fn precompute(op: &mut Op) {
-    let ellps = op.params.ellps[0];
+    let ellps = op.params.ellps(0);
     let n = ellps.third_flattening();
-    let lat_0 = op.params.lat[0];
-    let y_0 = op.params.y[0];
+    let lat_0 = op.params.lat(0).to_radians();
+    let y_0 = op.params.y(0);
 
     // The scaled spherical Earth radius - Qn in Engsager's implementation
-    let qs = op.params.k[0] * ellps.semimajor_axis() * ellps.normalized_meridian_arc_unit();
+    let qs = op.params.k(0) * ellps.semimajor_axis() * ellps.normalized_meridian_arc_unit();
     op.params.real.insert("scaled_radius", qs);
 
     // The Fourier series for the conformal latitude
@@ -303,51 +305,10 @@ pub fn new(parameters: &RawParameters, ctx: &dyn Context) -> Result<Op, Error> {
 #[cfg(test)]
 mod tests {
     use super::*;
+    use float_eq::assert_float_eq;
 
     #[test]
     fn tmerc() -> Result<(), Error> {
-        // Validation values from PROJ:
-        // echo 12 55 0 0 | cct -d18 +proj=utm +zone=32 | clip
-        #[rustfmt::skip]
-        let geo = [
-            Coor4D::geo( 55.,  12., 0., 0.),
-            Coor4D::geo(-55.,  12., 0., 0.),
-            Coor4D::geo( 55., -6., 0., 0.),
-            Coor4D::geo(-55., -6., 0., 0.)
-        ];
-
-        #[rustfmt::skip]
-        let projected = [
-            Coor4D::raw( 691_875.632_139_661, 6_098_907.825_005_012, 0., 0.),
-            Coor4D::raw( 691_875.632_139_661,-6_098_907.825_005_012, 0., 0.),
-            Coor4D::raw(-455_673.814_189_040, 6_198_246.671_090_279, 0., 0.),
-            Coor4D::raw(-455_673.814_189_040,-6_198_246.671_090_279, 0., 0.)
-        ];
-
-        let mut ctx = Minimal::default();
-        let definition = "tmerc k_0=0.9996 lon_0=9 x_0=500000";
-        let op = ctx.op(definition)?;
-
-        let mut operands = geo.clone();
-        ctx.apply(op, Fwd, &mut operands)?;
-
-        for i in 0..operands.len() {
-            dbg!(operands[i]);
-            dbg!(projected[i]);
-            assert!(operands[i].hypot2(&projected[i]) < 1e-6);
-        }
-
-        ctx.apply(op, Inv, &mut operands)?;
-        for i in 0..operands.len() {
-            assert!(operands[i].hypot2(&geo[i]) < 5e-6);
-        }
-
-        Ok(())
-    }
-
-    // Same as above, but using the 2D Coor2D data structure
-    #[test]
-    fn coor2d() -> Result<(), Error> {
         // Validation values from PROJ:
         // echo 12 55 0 0 | cct -d18 +proj=utm +zone=32 | clip
         #[rustfmt::skip]
@@ -374,9 +335,7 @@ mod tests {
         ctx.apply(op, Fwd, &mut operands)?;
 
         for i in 0..operands.len() {
-            dbg!(operands[i]);
-            dbg!(projected[i]);
-            assert!(operands[i].hypot2(&projected[i]) < 1e-6);
+            assert_float_eq!(operands[i].0, projected[i].0, abs_all <= 1e-8);
         }
 
         ctx.apply(op, Inv, &mut operands)?;
@@ -397,29 +356,29 @@ mod tests {
         // echo 12 55 0 0 | cct -d18 +proj=utm +zone=32 | clip
         #[rustfmt::skip]
         let geo = [
-            Coor4D::geo( 55.,  12., 0., 0.),
-            Coor4D::geo(-55.,  12., 0., 0.),
-            Coor4D::geo( 55., -6., 0., 0.),
-            Coor4D::geo(-55., -6., 0., 0.)
+            Coor2D::geo( 55.,  12.),
+            Coor2D::geo(-55.,  12.),
+            Coor2D::geo( 55., -6.),
+            Coor2D::geo(-55., -6.)
         ];
 
         #[rustfmt::skip]
         let projected = [
-            Coor4D::raw( 691_875.632_139_661, 6_098_907.825_005_012, 0., 0.),
-            Coor4D::raw( 691_875.632_139_661,-6_098_907.825_005_012, 0., 0.),
-            Coor4D::raw(-455_673.814_189_040, 6_198_246.671_090_279, 0., 0.),
-            Coor4D::raw(-455_673.814_189_040,-6_198_246.671_090_279, 0., 0.)
+            Coor2D::raw( 691_875.632_139_661, 6_098_907.825_005_012),
+            Coor2D::raw( 691_875.632_139_661,-6_098_907.825_005_012),
+            Coor2D::raw(-455_673.814_189_040, 6_198_246.671_090_279),
+            Coor2D::raw(-455_673.814_189_040,-6_198_246.671_090_279)
         ];
 
         let mut operands = geo.clone();
-        ctx.apply(op, Fwd, &mut operands)?;
+        assert_eq!(ctx.apply(op, Fwd, &mut operands)?, 4);
         for i in 0..operands.len() {
-            assert!(operands[i].hypot2(&projected[i]) < 5e-3);
+            assert_float_eq!(operands[i].0, projected[i].0, abs_all <= 1e-8);
         }
 
-        ctx.apply(op, Inv, &mut operands)?;
+        assert_eq!(ctx.apply(op, Inv, &mut operands)?, 4);
         for i in 0..operands.len() {
-            assert!(operands[i].hypot2(&geo[i]) < 10e-8);
+            assert_float_eq!(operands[i].0, geo[i].0, abs_all <= 1e-12);
         }
 
         Ok(())
diff --git a/src/inner_op/webmerc.rs b/src/inner_op/webmerc.rs
index cbe4186..0da3bde 100644
--- a/src/inner_op/webmerc.rs
+++ b/src/inner_op/webmerc.rs
@@ -6,7 +6,7 @@ use std::f64::consts::FRAC_PI_4;
 // ----- F O R W A R D -----------------------------------------------------------------
 
 fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
-    let ellps = op.params.ellps[0];
+    let ellps = op.params.ellps(0);
     let a = ellps.semimajor_axis();
 
     let mut successes = 0_usize;
@@ -30,7 +30,7 @@ fn fwd(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
 // ----- I N V E R S E -----------------------------------------------------------------
 
 fn inv(op: &Op, _ctx: &dyn Context, operands: &mut dyn CoordinateSet) -> usize {
-    let ellps = op.params.ellps[0];
+    let ellps = op.params.ellps(0);
     let a = ellps.semimajor_axis();
 
     let mut successes = 0_usize;
@@ -80,6 +80,7 @@ pub fn new(parameters: &RawParameters, _ctx: &dyn Context) -> Result<Op, Error>
 #[cfg(test)]
 mod tests {
     use super::*;
+    use float_eq::assert_float_eq;
 
     #[test]
     fn webmerc() -> Result<(), Error> {
@@ -102,14 +103,13 @@ mod tests {
         let mut operands = geo.clone();
         ctx.apply(op, Fwd, &mut operands)?;
         for i in 0..operands.len() {
-            assert!(operands[i].hypot2(&projected[i]) < 1e-8);
+            assert_float_eq!(operands[i].0, projected[i].0, abs_all <= 1e-8);
         }
 
         // Roundtrip
         ctx.apply(op, Inv, &mut operands)?;
         for i in 0..operands.len() {
-            dbg!(operands[0].to_geo());
-            assert!(operands[i].hypot2(&geo[i]) < 2e-9);
+            assert_float_eq!(operands[i].0, geo[i].0, abs_all <= 2e-9);
         }
 
         Ok(())
diff --git a/src/op/mod.rs b/src/op/mod.rs
index 4689216..0a62e72 100644
--- a/src/op/mod.rs
+++ b/src/op/mod.rs
@@ -67,7 +67,19 @@ impl Op {
         _ctx: &dyn Context,
     ) -> Result<Op, Error> {
         let def = parameters.definition.as_str();
-        let params = ParsedParameters::new(parameters, gamut)?;
+        let mut params = ParsedParameters::new(parameters, gamut)?;
+
+        // Convert lat_{0..4} and lon_{0..4} to radians
+        for i in ["lat_0", "lat_1", "lat_2", "lat_3"] {
+            let lat = *params.real.get(i).unwrap_or(&0.);
+            params.real.insert(i, lat);
+        }
+
+        for i in ["lon_0", "lon_1", "lon_2", "lon_3"] {
+            let lon = *params.real.get(i).unwrap_or(&0.);
+            params.real.insert(i, lon);
+        }
+
         let descriptor = OpDescriptor::new(def, fwd, inv);
         let steps = Vec::<Op>::new();
         let id = OpHandle::new();
diff --git a/src/op/parameter.rs b/src/op/parameter.rs
index 0344158..3ddc7dc 100644
--- a/src/op/parameter.rs
+++ b/src/op/parameter.rs
@@ -12,7 +12,7 @@
 ///
 /// For a given operation, the union of the sets of its required and optional
 /// parameters is called the *gamut* of the operation.
-#[derive(Debug)]
+#[derive(Debug, PartialEq, Clone)]
 pub enum OpParameter {
     /// A flag is a boolean that is true if present, false if not
     Flag { key: &'static str },
diff --git a/src/op/parsed_parameters.rs b/src/op/parsed_parameters.rs
index d7c57f8..b64daf4 100644
--- a/src/op/parsed_parameters.rs
+++ b/src/op/parsed_parameters.rs
@@ -11,18 +11,23 @@ use std::collections::BTreeSet;
 
 use super::*;
 
+#[rustfmt::skip]
+const ZERO_VALUED_IMPLICIT_GAMUT_ELEMENTS: [&str; 16] = [
+    "x_0", "x_1", "x_2", "x_3",
+    "y_0", "y_1", "y_2", "y_3",
+    "lat_0", "lat_1", "lat_2", "lat_3",
+    "lon_0", "lon_1", "lon_2", "lon_3"
+];
+
+#[rustfmt::skip]
+const UNIT_VALUED_IMPLICIT_GAMUT_ELEMENTS: [&str; 4] = [
+    "k_0", "k_1", "k_2", "k_3"
+];
+
 #[derive(Debug, Clone)]
 pub struct ParsedParameters {
     pub name: String,
 
-    // Commonly used options have hard-coded slots
-    pub ellps: [Ellipsoid; 2],
-    pub lat: [f64; 4],
-    pub lon: [f64; 4],
-    pub x: [f64; 4],
-    pub y: [f64; 4],
-    pub k: [f64; 4],
-
     // Op-specific options are stored in B-Trees
     pub boolean: BTreeSet<&'static str>,
     pub natural: BTreeMap<&'static str, usize>,
@@ -88,23 +93,34 @@ impl ParsedParameters {
     pub fn ignored(&self) -> Vec<String> {
         self.ignored.clone()
     }
-    pub fn ellps(&self, index: usize) -> &Ellipsoid {
-        &self.ellps[index]
+    pub fn ellps(&self, index: usize) -> Ellipsoid {
+        // if 'ellps' was explicitly given, it will override 'ellps_0'
+        if index == 0 {
+            if let Some(e) = self.text.get("ellps") {
+                return Ellipsoid::named(e).unwrap();
+            }
+        }
+        let key = format!("ellps_{index}");
+        if let Some(e) = self.text.get(&key[..]) {
+            return Ellipsoid::named(e).unwrap();
+        }
+        // If none of them existed, i.e. no defaults were given, we return the general default
+        Ellipsoid::default()
+    }
+    pub fn k(&self, index: usize) -> f64 {
+        *(self.real.get(&format!("k_{index}")[..]).unwrap_or(&1.))
     }
     pub fn x(&self, index: usize) -> f64 {
-        self.x[index]
+        *(self.real.get(&format!("x_{index}")[..]).unwrap_or(&0.))
     }
     pub fn y(&self, index: usize) -> f64 {
-        self.y[index]
+        *(self.real.get(&format!("y_{index}")[..]).unwrap_or(&0.))
     }
     pub fn lat(&self, index: usize) -> f64 {
-        self.lat[index]
+        *self.real.get(&format!("lat_{index}")[..]).unwrap_or(&0.)
     }
     pub fn lon(&self, index: usize) -> f64 {
-        self.lon[index]
-    }
-    pub fn k(&self, index: usize) -> f64 {
-        self.k[index]
+        *self.real.get(&format!("lon_{index}")[..]).unwrap_or(&0.)
     }
 }
 
@@ -127,18 +143,9 @@ impl ParsedParameters {
         let mut uuid = BTreeMap::<&'static str, uuid::Uuid>::new();
         let fourier_coefficients = BTreeMap::<&'static str, FourierCoefficients>::new();
 
-        // 'omit_fwd'/'omit_inv' are implicitly valid for all operators (including pipelines),
-        // so we append them to the gamut
-        let mut gamutt = Vec::new();
-        for g in gamut {
-            gamutt.push(g);
-        }
-        gamutt.push(&OpParameter::Flag { key: "omit_fwd" });
-        gamutt.push(&OpParameter::Flag { key: "omit_inv" });
-
         // Try to locate all accepted parameters, type check, and place them into
         // their proper bins
-        for p in gamutt {
+        for p in gamut {
             match *p {
                 OpParameter::Flag { key } => {
                     if let Some(value) = chase(globals, &locals, key)? {
@@ -311,55 +318,30 @@ impl ParsedParameters {
             };
         }
 
-        // Now handle the commonly used options with the hard-coded slots
+        // Default gamut elements - traditionally supported for all operators
 
-        let mut ellps = [Ellipsoid::default(), Ellipsoid::default()];
-        let mut lat = [0.; 4];
-        let mut lon = [0.; 4];
-        let mut x = [0.; 4];
-        let mut y = [0.; 4];
-        let mut k = [1.; 4];
-
-        // ellps_{n}
-        for i in 0..2 {
-            let key = format!("ellps_{i}");
-            if let Some(e) = text.get(&key[..]) {
-                ellps[i] = Ellipsoid::named(e)?;
+        // omit_fwd and omit_inv are implicitly valid for all ops
+        if let Some(value) = chase(globals, &locals, "omit_fwd")? {
+            if value.is_empty() || value.to_lowercase() == "true" {
+                boolean.insert("omit_fwd");
             }
         }
-        // But `ellps` trumps `ellps_0`
-        if let Some(e) = text.get("ellps") {
-            ellps[0] = Ellipsoid::named(e)?;
-        }
-
-        // lat_{n}
-        for i in 0..4 {
-            let key = format!("lat_{i}");
-            lat[i] = (*real.get(&key[..]).unwrap_or(&0.)).to_radians();
-        }
-
-        // lon_{n}
-        for i in 0..4 {
-            let key = format!("lon_{i}");
-            lon[i] = (*real.get(&key[..]).unwrap_or(&0.)).to_radians();
-        }
-
-        // x_{n}
-        for i in 0..4 {
-            let key = format!("x_{i}");
-            x[i] = *real.get(&key[..]).unwrap_or(&0.);
+        if let Some(value) = chase(globals, &locals, "omit_inv")? {
+            if value.is_empty() || value.to_lowercase() == "true" {
+                boolean.insert("omit_inv");
+            }
         }
 
-        // y_{n}
-        for i in 0..4 {
-            let key = format!("y_{i}");
-            y[i] = *real.get(&key[..]).unwrap_or(&0.);
+        for k in ZERO_VALUED_IMPLICIT_GAMUT_ELEMENTS {
+            if !real.contains_key(k) {
+                real.insert(k, 0.);
+            }
         }
 
-        // k_{n}
-        for i in 0..4 {
-            let key = format!("k_{i}");
-            k[i] = *real.get(&key[..]).unwrap_or(&0.);
+        for k in UNIT_VALUED_IMPLICIT_GAMUT_ELEMENTS {
+            if !real.contains_key(k) {
+                real.insert(k, 1.);
+            }
         }
 
         let name = locals
@@ -374,12 +356,6 @@ impl ParsedParameters {
         let given = locals.clone();
         let ignored: Vec<String> = locals.into_keys().collect();
         Ok(ParsedParameters {
-            ellps,
-            lat,
-            lon,
-            x,
-            y,
-            k,
             name,
             boolean,
             natural,
@@ -522,7 +498,7 @@ mod tests {
         assert_eq!(*p.text.get("text").unwrap(), "text");
 
         assert_eq!(
-            p.ellps[0].semimajor_axis(),
+            p.ellps(0).semimajor_axis(),
             Ellipsoid::new(123., 1. / 456.).semimajor_axis()
         );