use std::{cell::RefCell, collections::{BTreeMap, HashMap}, num::NonZero, ops::Deref, rc::Rc};
use colored::Colorize;

use anyhow::Result;
use kiddo::{Manhattan, NearestNeighbour, float::kdtree::KdTree};
use log::{debug, info, trace, warn};
use nalgebra::{Point3, Unit, Vector3};
use obj::ObjData;

use crate::{
    Aabb, Collide, IdxTriangle, Split, SplitEdges, primitives::{IdxEdge, IdxIntersection, Normal, PrimitiveIdx, Vertex}
};

type Dir3 = Unit<Vector3<f64>>;

#[derive(Clone, Debug)]
pub struct Verts {
    verts: Rc<RefCell<Vec<Point3<f64>>>>,
}

impl Deref for Verts {
    type Target = Rc<RefCell<Vec<Point3<f64>>>>;
    fn deref(&self) -> &Rc<RefCell<Vec<Point3<f64>>>> {
        &self.verts
    }
}

impl Verts {
    pub fn new() -> Self {
        let verts = Rc::new(RefCell::new(Vec::new()));
        Self { verts }
    }
    pub fn with_capacity(capacity: usize) -> Self {
        let verts = Rc::new(RefCell::new(Vec::with_capacity(capacity)));
        Self { verts }
    }
    pub fn alllocate_size(&self, size: usize) {
        let current_size = self.verts.borrow().len();
        self.verts.borrow_mut().resize(current_size + size, Point3::new(0.0, 0.0, 0.0));
    }
    pub fn allocate(&self, vert: Point3<f64>) -> PrimitiveIdx {
        let idx = self.verts.borrow().len();
        self.verts.borrow_mut().push(vert);
        PrimitiveIdx::Global(idx)
    }
    pub fn allocate_iter(&self, verts: impl IntoIterator<Item=Point3<f64>>) -> Vec<PrimitiveIdx> {
        let base_idx = self.verts.borrow().len();
        verts.into_iter().enumerate().map(|(i, vert)| {
            self.verts.borrow_mut().push(vert);
            PrimitiveIdx::Global(base_idx + i)
        }).collect()
    }
}

#[derive(Clone, Debug)]
pub struct Norms {
    norms: Rc<RefCell<Vec<Dir3>>>,
}

impl Deref for Norms {
    type Target = Rc<RefCell<Vec<Dir3>>>;
    fn deref(&self) -> &Rc<RefCell<Vec<Dir3>>> {
        &self.norms
    }
}

impl Norms {
    pub fn new() -> Self {
        let norms = Rc::new(RefCell::new(Vec::new()));
        Self { norms }
    }
    pub fn with_capacity(capacity: usize) -> Self {
        let norms = Rc::new(RefCell::new(Vec::with_capacity(capacity)));
        Self { norms }
    }
    pub fn alllocate_size(&self, size: usize) {
        let current_size = self.norms.borrow().len();
        self.norms.borrow_mut().resize(current_size + size, Unit::new_normalize(Vector3::new(1.0, 0.0, 0.0)));
    }
    pub fn allocate(&self, dir: Dir3) -> PrimitiveIdx {
        let idx = self.norms.borrow().len();
        self.norms.borrow_mut().push(dir);
        PrimitiveIdx::Global(idx)
    }
    pub fn allocate_iter(&self, dirs: impl IntoIterator<Item=Dir3>) -> Vec<PrimitiveIdx> {
        let base_idx = self.norms.borrow().len();
        dirs.into_iter().enumerate().map(|(i, dir)| {
            self.norms.borrow_mut().push(dir);
            PrimitiveIdx::Global(base_idx + i)
        }).collect()
    }
}

#[derive(Clone, Debug)]
pub struct Faces {
    faces: Rc<RefCell<Vec<IdxTriangle>>>,
}

impl Deref for Faces {
    type Target = Rc<RefCell<Vec<IdxTriangle>>>;
    fn deref(&self) -> &Rc<RefCell<Vec<IdxTriangle>>> {
        &self.faces
    }
}

impl Faces {
    pub fn new() -> Self {
        let faces = Rc::new(RefCell::new(Vec::new()));
        Self { faces }
    }
    pub fn with_capacity(capacity: usize) -> Self {
        let faces = Rc::new(RefCell::new(Vec::with_capacity(capacity)));
        Self { faces }
    }
    // FIXME: Allocate size should increase capacity, not resize and initialize to default
    //pub fn allocate_size(&self, size: usize) {
    //    let current_size = self.faces.borrow().len();
    //    self.faces.borrow_mut().resize(current_size + size, IdxTriangle::default());
    //}
    pub fn allocate(&self, tri: IdxTriangle) -> PrimitiveIdx {
        let idx = self.faces.borrow().len();
        self.faces.borrow_mut().push(tri);
        PrimitiveIdx::Global(idx)
    }
    pub fn allocate_iter(&self, tris: impl IntoIterator<Item=IdxTriangle>) -> Vec<PrimitiveIdx> {
        let base_idx = self.faces.borrow().len();
        tris.into_iter().enumerate().map(|(i, tri)| {
            self.faces.borrow_mut().push(tri);
            PrimitiveIdx::Global(base_idx + i)
        }).collect()
    }
    pub fn remap(&self, idx_rename: &BTreeMap<usize, usize>) {
        self.faces.borrow_mut().iter_mut().for_each(|tri| {
            tri.verts.iter_mut().for_each(|vert| {
                if let PrimitiveIdx::Global(idx) = vert.idx {
                    if let Some(new_idx) = idx_rename.get(&idx) {
                        vert.idx = PrimitiveIdx::Global(*new_idx);
                    }
                }
            });
        });
    }
}

pub fn parse_obj(
    obj_data: &ObjData,
) -> (Vec<Mesh>, Verts, Norms, Faces) {
    let verts = Verts::with_capacity(obj_data.position.len());
    let norms = Norms::with_capacity(obj_data.normal.len());
    let faces = Faces::new();

    let allocated_verts_idx = verts.allocate_iter(
        obj_data
        .position
        .iter()
        .map(|vs| {
            let vs_f64: [f64; 3] = vs.map(|v| v as f64);
            Point3::new(vs_f64[0], vs_f64[1], vs_f64[2])
        })
    );
    debug!("Allocated verts idx: {:?}", allocated_verts_idx);

    let allocated_norms_idx = norms.allocate_iter(
        obj_data
        .normal
        .iter()
        .map(|vs| {
            let vs_f64: [f64; 3] = vs.map(|v| v as f64);
            Dir3::new_normalize(Vector3::new(vs_f64[0], vs_f64[1], vs_f64[2]))
        })
    );

    debug!("Allocated norms idx: {:?}", allocated_norms_idx);

    let meshes =
        obj_data
        .objects
        .iter()
        .enumerate()
        .map(|(idx, obj)| {
            // NOTE: Collapse groups from an object
            let obj_faces: Vec<_> = obj
                .groups
                .iter()
                .flat_map(|group| {
                    group.polys.iter().map(|poly| {
                        let tri_verts: [Vertex; 3] = poly
                            .0
                            .iter()
                            .map(|idx_tuple| idx_tuple.0)
                            .map(|idx| Vertex::new(idx, &verts.borrow()))
                            .collect::<Vec<_>>()
                            .try_into()
                            .unwrap();
                        //let t_idx = poly.0.map(|idx_tuple| idx_tuple.1);
                        let tri_norms: Option<[Normal; 3]> = poly
                            .0
                            .iter()
                            .map(|idx_tuple| idx_tuple.2)
                            .map(|idx| idx.map(|i| Normal::new(i, &norms.borrow())))
                            .collect::<Option<Vec<_>>>()
                            .and_then(|v| v.try_into().ok());
                        (tri_verts, tri_norms)
                    })
                })
                .collect();

            // Tag IdxTriangle with the correct idx that would be allocated to faces
            let base_idx = faces.borrow().len();
            let idx_tris: Vec<_> = obj_faces
                .into_iter()
                .enumerate()
                .map(|(i, (tri_verts, tri_norms))| IdxTriangle::new(tri_verts, tri_norms, PrimitiveIdx::Global(base_idx + i)))
                .collect();
            //faces.alllocate_size(obj_faces.len());
            let polygons = faces.allocate_iter(idx_tris);
            debug!("Allocated faces idx for mesh {}: {:?}", obj.name, polygons);

            Mesh::from_polygons(
                obj.name.clone(),
                PrimitiveIdx::Global(idx),
                polygons,
                &verts,
                &norms,
                &faces)
        })
        .collect();

    // Vertex deduplication, and remap the vertex idx in faces accordingly
    let verts_remap = prune_verts(&verts, &faces);
    info!("Vertices remaped: {:?}", verts_remap.len());
    trace!("Vertices remaped: {:?}", verts_remap);

    (meshes, verts, norms, faces)
}

pub fn prune_verts(verts: &Verts, faces: &Faces) -> BTreeMap<usize, usize> {
    // Build k-d tree from vertices
    let mut kdtree: KdTree<f64, usize, 3, 128, u32> = KdTree::new();
    for (idx, vert) in verts.borrow().iter().enumerate() {
        kdtree.add(&[vert.x, vert.y, vert.z], idx);
    }
    // Map PrimitiveIdx::Global(idx) -> new PrimitiveIdx::Global(idx) after vertex deduplication
    let mut vert_rename: BTreeMap<usize, usize> = BTreeMap::new();
    for (idx, vert) in verts.borrow().iter().enumerate() {
        if vert_rename.contains_key(&idx) {
            continue;
        }
        // FIXME: Make the distance threshold for deduplication proportional to the size of the
        // scene
        let search_result = kdtree.nearest_n_within::<Manhattan>(&[vert.x, vert.y, vert.z], 1e-9, NonZero::new(10).unwrap(), true);
        for &NearestNeighbour{item: dup_idx, ..} in search_result.iter() {
            if dup_idx != idx {
                vert_rename.insert(dup_idx, idx);
            }
        }
    }
    faces.remap(&vert_rename);

    vert_rename
}

#[derive(Clone, Debug)]
pub struct Mesh {
    /// List of indexes pointing to the polygons in the list of shapes
    pub name: String,
    pub idx: PrimitiveIdx,
    pub polygons: Vec<PrimitiveIdx>,
    pub verts: Verts,
    pub norms: Norms,
    pub faces: Faces,
    pub aabb: Aabb,
}

impl Mesh {
    pub fn new(name: String, idx: PrimitiveIdx, verts: &Verts, norms: &Norms, faces: &Faces) -> Self {
        let aabb = Aabb::null();
        Self { name, idx, polygons: Vec::new(), verts: verts.clone(), norms: norms.clone(), faces: faces.clone(), aabb }
    }
    pub fn from_polygons(name: String, idx: PrimitiveIdx, polygons: Vec<PrimitiveIdx>, verts: &Verts, norms: &Norms, faces: &Faces) -> Self {
        let mut mesh = Self::new(name, idx, verts, norms, faces);
        let aabb = polygons.iter().fold(Aabb::null(), |aabb, &idx| {
            let tri_aabb = faces.borrow()[match idx {
                PrimitiveIdx::Local(_) => panic!("Expected global indices for mesh polygons, found: {:?}", idx),
                PrimitiveIdx::Global(idx) => idx,
            }].aabb();
            aabb + tri_aabb
        });
        mesh.polygons = polygons;
        mesh.aabb = aabb;
        mesh
    }
    pub fn push(&mut self, new_tris: Vec<IdxTriangle>) {
        let current_size = self.faces.borrow().len();
        for tri in new_tris.iter() {
            self.aabb = self.aabb.clone() + tri.aabb();
        }
        let new_polygons = self.faces.allocate_iter(
            new_tris
                .into_iter()
                .enumerate()
                .map(|(i, tri)| tri.into_global(current_size + i))
        );
        self.polygons.extend(new_polygons);
    }
    pub fn tris(&self) -> Vec<IdxTriangle> {
        self.polygons
            .iter()
            .map(|&idx| {
                let idx = match idx {
                    PrimitiveIdx::Local(_) => panic!("Expected global indices for mesh polygons, found: {:?}", idx),
                    PrimitiveIdx::Global(idx) => idx,
                };
                self.faces.borrow()[idx].clone()
            })
            .collect()
    }
    pub fn invert_normals(&self) {
        for tri in self.tris().iter_mut() {
            tri.invert_normals();
        }
    }
}

#[derive(Clone, Debug)]
pub struct FaceIntersection {
    pub tri_u: PrimitiveIdx,
    pub tri_v: PrimitiveIdx,
    pub split_edges: SplitEdges,
}

impl FaceIntersection {
    pub fn new(tri_u_idx: PrimitiveIdx, tri_v_idx: PrimitiveIdx, split_edges: SplitEdges) -> Self {
        let min_tri_idx = tri_u_idx.min(tri_v_idx);
        let max_tri_idx = tri_u_idx.max(tri_v_idx);
        assert!(min_tri_idx != max_tri_idx,
            "Intersection should be between two different triangles, but got the same triangle idx {:?}",
            tri_u_idx);
        Self { tri_u: min_tri_idx, tri_v: max_tri_idx, split_edges }
    }
    pub fn rename_vertices(&mut self, idx_alloc_map: &HashMap<IdxIntersection, PrimitiveIdx>) {
        self.split_edges.rename_vertices(idx_alloc_map);
    }
    pub fn into_tris(self, faces: &Faces) -> Result<Vec<IdxTriangle>> {
        let tri_u_idx = match self.tri_u {
            PrimitiveIdx::Local(_) => panic!("Expected global indices for face intersection triangles, found: {:?}", self.tri_u),
            PrimitiveIdx::Global(idx) => idx,
        };
        let tri_u = &faces.borrow()[tri_u_idx];
        tri_u.from_edges(self.split_edges.edges)
    }
}

impl PartialEq for FaceIntersection {
    fn eq(&self, other: &Self) -> bool {
        self.tri_u == other.tri_u && self.tri_v == other.tri_v
    }
}
impl Eq for FaceIntersection {}

impl Ord for FaceIntersection {
    fn cmp(&self, other: &Self) -> std::cmp::Ordering {
        match self.tri_u.cmp(&other.tri_u) {
            std::cmp::Ordering::Equal => self.tri_v.cmp(&other.tri_v),
            ord => ord,
        }
    }
}

impl PartialOrd for FaceIntersection {
    fn partial_cmp(&self, other: &Self) -> Option<std::cmp::Ordering> {
        Some(self.cmp(other))
    }
}

pub struct MeshIntersection {
    pub mesh_u: PrimitiveIdx,
    pub mesh_v: PrimitiveIdx,
    pub tri_intersections: Vec<FaceIntersection>,
}

pub struct VertexMatch {
    pub meshes: (PrimitiveIdx, PrimitiveIdx),
    pub tris: (PrimitiveIdx, PrimitiveIdx),
    pub from: (IdxEdge, IdxEdge),
}

// FIXME: A better colision detection could be employed to ensure
// the expensive Split methods are not called unecessary
impl Collide<Mesh> for Mesh {
    fn overlap(&self, other: &Mesh) -> bool {
        if !self.aabb.overlap(&other.aabb) {
            return false;
        }
        let mut tri_u_found = false;
        for &tri_u_idx in self.polygons.iter() {
            let tri_u = &self.faces.borrow()[*tri_u_idx];
            if tri_u.tri().overlap(&other.aabb) {
                tri_u_found = true;
                break;
            }
        }
        let mut tri_v_found = false;
        for &tri_v_idx in other.polygons.iter() {
            let tri_v = &other.faces.borrow()[*tri_v_idx];
            if tri_v.tri().overlap(&self.aabb) {
                tri_v_found = true;
                break;
            }
        }
        tri_u_found && tri_v_found
    }
}

// TODO: Replace asserts with Err
impl Split<Mesh, Vec<FaceIntersection>> for Mesh {
    type Inst = Mesh;
    fn intersect(&self, other: &Mesh) -> Result<Vec<FaceIntersection>> {
        // If there is no intersection, no point to proceed with checking pairwise triangles
        if !self.overlap(other) {
            return Ok(Vec::new());
        }

        //  Extract the triangles from mesh, keeping track of their idx for regrouping
        //  intersections
        let tri_u_union: Vec<_> = self
            .polygons
            .iter()
            .map(|&idx| {
                (
                    idx,
                    self.faces.borrow()[*idx].clone(),
                )
            })
            .collect();
        let tri_v_union: Vec<_> = other
            .polygons
            .iter()
            .map(|&idx| {
                (
                    idx,
                    other.faces.borrow()[*idx].clone(),
                )
            })
            .collect();

        // Find all intersections
        let mut intersections: Vec<FaceIntersection> = Vec::new();
        for (idx_tri_u, tri_u) in tri_u_union.iter() {
            for (idx_tri_v, tri_v) in tri_v_union.iter() {
                // FIXME: Would be better to have a collision tree checking here instead of
                // brute force pairwise collision tests
                if tri_u.overlap(tri_v) {
                    debug!("Intersect triangles {:?} and {:?} from meshes {:?} and {:?}",
                        idx_tri_u, idx_tri_v, self.idx, other.idx);
                    let new_intersection = tri_u.intersect(tri_v)?;
                    if new_intersection.edges.is_empty() {
                        warn!("Triangles {:?} and {:?} from meshes {:?} and {:?} overlap but do not intersect",
                            idx_tri_u, idx_tri_v, self.idx, other.idx);
                    } else {
                        intersections.push(FaceIntersection::new(
                            *idx_tri_u,
                            *idx_tri_v,
                            new_intersection,
                        ));
                    }
                }
            }
        }
        info!("Found {} intersecting triangle pairs between meshes {:?} and {:?}",
            intersections.len(), self.idx, other.idx);
        Ok(intersections)
    }

    fn split(&self, other: Vec<FaceIntersection>) -> Result<(Mesh, Vec<FaceIntersection>)> {
        let tri_u_union: Vec<_> = self
            .polygons
            .iter()
            .map(|&idx| {
                (
                    idx,
                    self.faces.borrow()[*idx].clone(),
                )
            })
            .collect();

        let mut new_mesh = Mesh {
            name: self.name.clone(),
            idx: self.idx,
            polygons: Vec::new(),
            verts: self.verts.clone(),
            norms: self.norms.clone(),
            faces: self.faces.clone(),
            aabb: Aabb::null(),
        };

        // Map from local vertex idx in intersection to new global vertex idx in new_mesh
        let mut vert_idx_alloc_map: HashMap<IdxIntersection, PrimitiveIdx> = HashMap::new();
        let mut norm_idx_alloc_map: HashMap<PrimitiveIdx, PrimitiveIdx> = HashMap::new();

        // Collect split edges derived from the target triangle
        for (idx_tri_u, tri_u) in tri_u_union.iter() {
            // Collect all splits with tri_u, identified by idx_tri_u:PrimitiveIdx
            let mut splits = SplitEdges::new(vec![]);
            // TODO: Instead of iterating through all, sort intersection and binary search for the
            // ones matching idx_tri_u
            other.iter().for_each(|intersection| {
                if intersection.tri_u == *idx_tri_u || intersection.tri_v == *idx_tri_u {
                    splits = splits.union(&intersection.split_edges);
                }
            });

            // Allocate vertices from splits, and keep track of the mapping from local vertex idx to global
            let local_verts = Vertex::from_edges(&splits.edges)
                .into_iter()
                .filter(|vert| matches!(vert.idx, PrimitiveIdx::Local(..)));

            for vert in local_verts {
                assert!(vert.from.is_some(), "Local vertex must have intersection edges described");
                if !vert_idx_alloc_map.contains_key(&vert.from.unwrap()) {
                    let new_idx = new_mesh.verts.allocate(vert.value);
                    debug_assert!(matches!(new_idx, PrimitiveIdx::Global(_)), "Vertices from edges should be allocated globally");
                    vert_idx_alloc_map.insert(vert.from.unwrap(), new_idx);
                }
            }

            // Rneame vertices idx, mapping to the allocated global idx
            splits.rename_vertices(&vert_idx_alloc_map);

            // Sanity check: All vertices in splits should now be allocated globally
            Vertex::from_edges(&splits.edges)
                .into_iter()
                .for_each(|vert| {
                    assert!(matches!(vert.idx, PrimitiveIdx::Global(_)), "All vertices in splits should be allocated globally");
                });

            if splits.is_empty() {
                // No splits, push the original triangle to the new mesh
                new_mesh.push(vec![tri_u.clone()]);
            } else {
                // Split the triangle with the splits, and push to the new mesh
                let (mut new_tris_u, _splits) = tri_u.split(splits)?;

                // Rename norms in new_tris_u, mapping to the allocated global idx
                // WARN: Normals are likely to be duplicates for curved surface, but they should only affect the exported format
                for tri in new_tris_u.iter_mut() {
                    if let Some(norms) = tri.norms {
                        for norm in norms {
                            if matches!(norm.idx, PrimitiveIdx::Local(_)) && !norm_idx_alloc_map.contains_key(&norm.idx) {
                                let new_idx = new_mesh.norms.allocate(norm.value);
                                debug_assert!(matches!(new_idx, PrimitiveIdx::Global(_)), "Normals from original triangle should be allocated globally");
                                norm_idx_alloc_map.insert(norm.idx, new_idx);
                            }
                        }
                    }
                    tri.rename_norms(&norm_idx_alloc_map);
                }

                new_mesh.push(new_tris_u);
            }
        }

        let new_face_inters = other
            .iter()
            .map(|face_inter| {
                let mut new_face_inter = face_inter.clone();
                new_face_inter.rename_vertices(&vert_idx_alloc_map);
                new_face_inter
            })
            .collect();

        Ok((new_mesh, new_face_inters))
    }
}

pub fn remesh(mut meshes: Vec<Mesh>) -> Result<Vec<Mesh>> {
    info!("Start splitting meshes");
    for mesh in &meshes {
        debug!(" > Mesh: {}", mesh.name.blue());
    }

    let mut resolved_meshes = Vec::new();
    let mut idx = meshes.len();
    for i in 0..meshes.len() {
        for j in i+1..meshes.len() {
            let mesh_a = &meshes[i];
            let mesh_b = &meshes[j];
            if mesh_a.overlap(mesh_b) {
                info!("Mesh {} and {} overlap", mesh_a.name, mesh_b.name);
            } else {
                continue;
            }
            let inter = mesh_a.intersect(mesh_b)?;
            trace!("Mesh {} and {} intersection: {:?}", mesh_a.name, mesh_b.name, inter);
            let (new_mesh_a, inter) = mesh_a.split(inter)?;
            let (new_mesh_b, inter) = mesh_b.split(inter)?;
            if !inter.is_empty() {
                let mut inter_mesh = Mesh {
                    name: format!("{}-{}", mesh_a.name, mesh_b.name),
                    idx: PrimitiveIdx::Global(idx),
                    polygons: Vec::new(),
                    verts: mesh_a.verts.clone(),
                    norms: mesh_a.norms.clone(),
                    faces: mesh_a.faces.clone(),
                    aabb: Aabb::null(),
                };
                for tri_inter in inter.into_iter() {
                    // WARN: Information is lost about the normals in the intersection, need to
                    // propagate them from the source triangles
                    inter_mesh.push(tri_inter.into_tris(&mesh_a.faces)?);
                };
                resolved_meshes.push(inter_mesh);
                idx += 1;
            }
            meshes[i] = new_mesh_a;
            meshes[j] = new_mesh_b;
        }
        resolved_meshes.push(meshes[i].clone());
    }

    Ok(resolved_meshes)
}

#[cfg(test)]
mod tests {
    use std::{io::BufReader, path::Path, sync::Once};
    use obj::ObjData;
    use tempfile::NamedTempFile;
    use super::*;
    use crate::save::Save;

    static INIT: Once = Once::new();
    fn init_logger() {
        INIT.call_once(|| {
            env_logger::Builder::from_default_env()
                .is_test(true) // don't add timestamps
                .try_init()
                .ok();
        });
    }

    #[test]
    fn test_parse_obj() -> Result<()> {
        init_logger();
        static OBJ_STR: &'static str = "
o Tank
v 0 0 0
v 2 0 0
v 2 1 0
v 0 1 0
vn 0 0 1
f 1//1 2//1 3//1
f 1//1 3//1 4//1

o Water
v 2 0.8 0
v 0 0.8 0
vn 0 0 -1
f 1//2 2//2 6//2
f 2//2 5//2 6//2
";
        let mut reader = BufReader::new(OBJ_STR.as_bytes());
        let obj = ObjData::load_buf(&mut reader)?;
        let (meshes, verts, norms, faces) = parse_obj(&obj);

        println!("Verts: {}", verts.borrow().len());
        println!("Norms: {}", norms.borrow().len());
        println!("Faces: {}", faces.borrow().len());

        let resolved_meshes = remesh(meshes)?;

        println!("Resolved Meshes: {:?}", resolved_meshes.iter().map(|m| m.polygons.clone()));

        assert!(resolved_meshes.len() == 3, "Expected 3 meshes after remeshing, got {}", resolved_meshes.len());

        assert_eq!(resolved_meshes[0].polygons.len(), 5, "Intersection mesh should have 5 triangles");
        assert_eq!(resolved_meshes[1].polygons.len(), 3, "Tank remained mesh should have 3 triangles");
        assert_eq!(resolved_meshes[2].polygons.len(), 0, "Water surface is contained within the tank surface");


        // Rexport the meshes to buffer and print to screen
        let infile = NamedTempFile::new().expect("Expected Temporary file to write test spectrum");
        let tmp_path_str = infile.path().to_str().unwrap().to_string();
        println!("Re-export the mesh to {}", tmp_path_str);

        resolved_meshes.save_data(Path::new(&tmp_path_str))?;

        // Read the re-exported file and print to screen
        //let re_exported_str = std::fs::read_to_string(infile.path())?;
        //println!("Re-exported OBJ:\n{}", re_exported_str);

        Ok(())
    }

}
