//-----------------------------------------------------------------------
// <copyright file="PLLinesCreatorDupuitTests.cs" company="Deltares">
//  Copyright (c) 2011 Deltares. All rights reserved.
// </copyright>
// <author>B.S.T.I.M. The</author>
// <email>tom.the@deltares.nl</email>
// <date>31-08-2011</date>
// <summary>Write and read files to interact with DikeFlow.exe</summary>
//-----------------------------------------------------------------------

using Deltares.Geometry;
using Deltares.Geotechnics.Soils;
using Deltares.Geotechnics.SurfaceLines;
using Deltares.Mathematics;

namespace Deltares.Dam.Data
{
    using System;
    using System.Collections.Generic;
    using System.Globalization;
    using System.IO;
    using System.Linq;
    using System.Text;
    using System.Text.RegularExpressions;
    using System.Threading;

    /// <summary>
    /// Exception class for DupuitCalculatorFileIO
    /// </summary>
    public class DupuitCalculatorFileIOException : ApplicationException
    {
        public DupuitCalculatorFileIOException(string message)
            : base(message)
        {
        }
    }

    /// <summary>
    /// Write and read files to interact with DikeFlow.exe
    /// </summary>
    public class DupuitCalculatorFileIO
    {
        public const double tolerance = 0.00001;
        public const double cFactor_mPerSec_2_mPerDay = 86400.0;
        public const double cFactor_mPerDay_2_mPerSec = 1 / cFactor_mPerSec_2_mPerDay;
        public const string soilLayersFilenameExt = ".A";
        public const string waterLevelFilenameExt = ".B";
        public const string precipitationFilenameExt = ".E";
        public const string materialFilenameExt = ".G";
        public const string monitorFilenameExt = ".H";
        public const string safetyFactorFilenameExt = ".J";
        public const string sectionsFilenameExt = ".K";
        public const string calculateFilenameExt = ".P";
        public const string piezoLinesFilesnameExt = ".C";
        public const string polderLevelFilenameExt = ".U";

        public const string layerTypeTop = "top";
        public const string layerTypeDike = "dike";
        public const string layerTypeCover = "cover";
        public const string layerTypeAquifer1 = "aquifer1";
        public const string layerTypeAquifer2 = "aquifer2";
        public const string layerTypeAquitard = "aquitard";

        public const string soillayersId = "soillayers";
        public const string waterlevelId = "waterlevel";
        public const string precipitationId = "precipitation";
        public const string materialId = "material";
        public const string safetyFactorId = "safety";
        public const string monitorId = "monitor";
        public const string sectionsId = "sections";
        public const string calculateId = "calculate";
        public const string piezolinesId = "piezolines";
        public const string polderlevelId = "polderwaterlevel";

        private Geometry2DData geometry2DData = null;
        private Geometry2DData adaptedGeometry2DData;
        private TimeSerie waterLevelTimeserie = null;
        private TimeSerie safetyFactorTimeserie = null;
        private TimeSerie polderLevelTimeserie = null;
        private SoilList soilList = null;
        private SurfaceLine2 surfaceLine = null;
        private DupuitCalculationOptions dupuitCalculationOptions = null;
        private DupuitPLLinesTimeSerie dupuitPLLinesTimeSerie = null;
        private bool isReverseLayerOrder = false;

        /// <summary>
        /// constructor
        /// </summary>
        public DupuitCalculatorFileIO()
        {

        }

        /// <summary>
        /// Throwhelper
        /// </summary>
        private void ThrowIfFalse(bool condition, string message)
        {
            if (!condition)
            {
                throw new DupuitCalculatorFileIOException(message);
            }
        }

        /// <summary>
        /// Property Geometry2DData
        /// </summary>
        public Geometry2DData Geometry2DData
        {
            get { return geometry2DData; }
            set { geometry2DData = value; }
        }

        /// <summary>
        /// Property WaterLevelTimeserie
        /// </summary>
        public TimeSerie WaterLevelTimeserie
        {
            get { return waterLevelTimeserie; }
            set { waterLevelTimeserie = value; }
        }

        /// <summary>
        /// Property PolderLevelTimeserie
        /// </summary>
        public TimeSerie PolderLevelTimeserie
        {
            get { return polderLevelTimeserie; }
            set { polderLevelTimeserie = value; }
        }

        /// <summary>
        /// Property SoilList
        /// </summary>
        public SoilList SoilList
        {
            get { return soilList; }
            set { soilList = value; }
        }

        /// <summary>
        /// Property SurfaceLine
        /// </summary>
        public SurfaceLine2 SurfaceLine
        {
            get { return surfaceLine; }
            set { surfaceLine = value; }
        }

        /// <summary>
        /// Property DupuitCalculationOptions
        /// </summary>
        public DupuitCalculationOptions DupuitCalculationOptions
        {
            get { return dupuitCalculationOptions; }
            set { dupuitCalculationOptions = value; }
        }

        /// <summary>
        /// Property DupuitPlLinesTimeSerie
        /// </summary>
        public DupuitPLLinesTimeSerie DupuitPlLinesTimeSerie
        {
            get { return dupuitPLLinesTimeSerie; }
            set { dupuitPLLinesTimeSerie = value; }
        }

        /// <summary>
        /// Property IsReverseLayerOrder
        /// Used to write the layers in reversed order
        /// Unfortunately D-Geo Stability and DikeFlow use different order of layers
        /// </summary>
        public bool IsReverseLayerOrder
        {
            get { return isReverseLayerOrder; }
            set { isReverseLayerOrder = value; }
        }

        /// <summary>
        /// Property SafetyFactorTimeserie
        /// </summary>
        public TimeSerie SafetyFactorTimeserie
        {
            get { return safetyFactorTimeserie; }
            set { safetyFactorTimeserie = value; }
        }

        /// <summary>
        /// Determine if boundary is part of dike
        /// The criterium is that the line should be both above toe of dike  inside and outside
        /// </summary>
        /// <param name="boundaryLine"></param>
        /// <returns></returns>
        private bool IsBoundaryPartOfDike(Geometry2DBoundaryLine boundaryLine)
        {
            var dikeToeAtRiver = surfaceLine.CharacteristicPoints.GetGeometryPoint(CharacteristicPointType.DikeToeAtRiver);
            double dikeOutsideX = dikeToeAtRiver.X;
            double dikeOutsideY = dikeToeAtRiver.Z;
            var dikeToeAtPolder = surfaceLine.CharacteristicPoints.GetGeometryPoint(CharacteristicPointType.DikeToeAtPolder);
            double dikeInsideX = dikeToeAtPolder.X;
            double dikeInsideY = dikeToeAtPolder.Z;
            double y1 = boundaryLine.ZFromX(dikeOutsideX);
            double y2 = boundaryLine.ZFromX(dikeInsideX);
            return ((y1 >= dikeOutsideY) && (y2 >= dikeInsideY));
        }

        /// <summary>
        /// Determine what the types of layer for DikeFlow
        /// </summary>
        /// <param name="geometry2DData"></param>
        private void AssignLayerTypes(ref Geometry2DData geometry2DData)
        {
            geometry2DData.SetLayertype(0, LayerType.Top);
            for (int layerIndex = 0; layerIndex < geometry2DData.LayerCount; layerIndex++)
            {
                Geometry2DLayer layer = geometry2DData.GetLayer(layerIndex);
                if (layerIndex == geometry2DData.LayerCount - 1)
                {
                    geometry2DData.SetLayertype(layerIndex, LayerType.Aquifer1);
                }
                else
                {
                    if (IsBoundaryPartOfDike(layer.boundaryLine))
                    {
                        geometry2DData.SetLayertype(layerIndex + 1, LayerType.Dike);
                    }
                    else
                    {
                        geometry2DData.SetLayertype(layerIndex + 1, LayerType.Cover);
                    }
                }
            }
            // Make sure there is always at least one layer with dike material
            geometry2DData.SetLayertype(1, LayerType.Dike);
        }

        /// <summary>
        /// Write all relevant inputfiles
        /// </summary>
        /// <param name="filenameBase"></param>
        public void WriteInputFiles(string filenameBase)
        {
            Thread.CurrentThread.CurrentCulture = CultureInfo.InvariantCulture;
            adaptedGeometry2DData = ReverseGeometry2DData(geometry2DData);
            AssignLayerTypes(ref adaptedGeometry2DData);
            WriteGeometry2DData(filenameBase);
            WriteWaterLevels(filenameBase);
            WritePolderLevels(filenameBase);
            WritePrecipitation(filenameBase);
            WriteMonitor(filenameBase);
            WriteMaterials(filenameBase);
            WriteCalculationOptions(filenameBase);
            WriteSections(filenameBase);
        }

        /// <summary>
        /// Read all relevant inputfiles (only needed for testing
        /// </summary>
        /// <param name="filenameBase"></param>
        public void ReadInputFiles(string filenameBase)
        {
            throw new NotImplementedException();
        }

        /// <summary>
        /// Write outputfiles (for testing reading outputfiles
        /// </summary>
        /// <param name="filenameBase"></param>
        public void WriteOutputFiles(string filenameBase)
        {
            DupuitPLLinesTimeSerie originalDupuitPLLinesTimeSerie;
            originalDupuitPLLinesTimeSerie = dupuitPLLinesTimeSerie;
            dupuitPLLinesTimeSerie = ReverseDupuitPLLinesTimeSerie(originalDupuitPLLinesTimeSerie);
            WritePlLines(filenameBase);
            dupuitPLLinesTimeSerie = originalDupuitPLLinesTimeSerie;

            WriteSafetyFactors(filenameBase);
        }

        /// <summary>
        /// Read relevant outputfiles
        /// </summary>
        /// <param name="filenameBase"></param>
        public void ReadOutputFiles(string filenameBase)
        {
            ReadPlLines(filenameBase);
            dupuitPLLinesTimeSerie = ReverseDupuitPLLinesTimeSerie(dupuitPLLinesTimeSerie);

            ReadSafetyFactors(filenameBase);
        }

        /// <summary>
        /// Write geometry data to file
        /// </summary>
        /// <param name="filenameBase"></param>
        private void WriteGeometry2DData(string filenameBase)
        {
            if (adaptedGeometry2DData != null)
            {
                var contentStringBuilder = new StringBuilder();
                contentStringBuilder.AppendLine(soillayersId);
                contentStringBuilder.AppendLine(String.Format("nline   {0}", adaptedGeometry2DData.LayerCount));
                for (int lineIndex = 0; lineIndex < adaptedGeometry2DData.LayerCount; lineIndex++)
                {
                    contentStringBuilder.AppendLine(String.Format("{0}  {1}  {2}",
                         adaptedGeometry2DData.GetLayer(lineIndex).boundaryLine.Points.Count, 2, lineIndex + 1));
                    for (int pointIndex = 0; pointIndex < adaptedGeometry2DData.GetLayer(lineIndex).boundaryLine.Points.Count; pointIndex++)
                    {
                        contentStringBuilder.AppendLine(String.Format("{0:e}  {1:e}",
                             this.adaptedGeometry2DData.GetLayer(lineIndex).boundaryLine.Points[pointIndex].X,
                             this.adaptedGeometry2DData.GetLayer(lineIndex).boundaryLine.Points[pointIndex].Z));
                    }
                }

                string contentString = contentStringBuilder.ToString();
                System.IO.File.WriteAllText(filenameBase + DupuitCalculatorFileIO.soilLayersFilenameExt, contentString);
            }

        }

        /// <summary>
        /// Write waterlevels to file
        /// </summary>
        /// <param name="filenameBase"></param>
        public void WriteWaterLevels(string filenameBase)
        {
            if (WaterLevelTimeserie != null)
            {
                var fullcontentStringBuilder = new StringBuilder();
                var contentStringBuilder = new StringBuilder();
                int currentTotal = 0;
                for (int waterlevelIndex = 0;
                     waterlevelIndex < this.WaterLevelTimeserie.Entries.Count();
                     waterlevelIndex++)
                {
                    if (!ComputeDoubles.IsNearEqual(this.WaterLevelTimeserie.Entries[waterlevelIndex].Value, this.WaterLevelTimeserie.MissVal, tolerance))
                    {
                        currentTotal++;
                        contentStringBuilder.AppendLine(String.Format("{0:e}  {1:e}",
                            (this.WaterLevelTimeserie.Entries[waterlevelIndex].DateTime - this.WaterLevelTimeserie.StartDateTime).TotalDays,
                            this.WaterLevelTimeserie.Entries[waterlevelIndex].Value));
                    }
                }
                fullcontentStringBuilder.AppendLine(waterlevelId);
                fullcontentStringBuilder.AppendLine(String.Format("{0}  {1}  {2}", currentTotal, 2, 1));
                fullcontentStringBuilder.Append(contentStringBuilder);

                string contentString = fullcontentStringBuilder.ToString();
                System.IO.File.WriteAllText(filenameBase + DupuitCalculatorFileIO.waterLevelFilenameExt, contentString);
            }
        }

        /// <summary>
        /// Write polder waterlevels to file
        /// </summary>
        /// <param name="filenameBase"></param>
        private void WritePolderLevels(string filenameBase)
        {
            if (PolderLevelTimeserie != null)
            {
                var contentStringBuilder = new StringBuilder();
                contentStringBuilder.AppendLine(polderlevelId);
                contentStringBuilder.AppendLine(String.Format("{0}  {1}  {2}", this.PolderLevelTimeserie.Entries.Count(), 2, 1));
                for (int polderlevelIndex = 0;
                     polderlevelIndex < this.PolderLevelTimeserie.Entries.Count();
                     polderlevelIndex++)
                {
                    contentStringBuilder.AppendLine(String.Format("{0:e}  {1:e}",
                        (this.PolderLevelTimeserie.Entries[polderlevelIndex].DateTime - this.PolderLevelTimeserie.StartDateTime).TotalDays,
                        this.polderLevelTimeserie.Entries[polderlevelIndex].Value));
                }

                string contentString = contentStringBuilder.ToString();
                System.IO.File.WriteAllText(filenameBase + DupuitCalculatorFileIO.polderLevelFilenameExt, contentString);
            }
        }

        /// <summary>
        /// Write safety to file
        /// </summary>
        /// <param name="filenameBase"></param>
        private void WriteSafetyFactors(string filenameBase)
        {
            if (WaterLevelTimeserie != null)
            {
                var contentStringBuilder = new StringBuilder();
                contentStringBuilder.AppendLine(safetyFactorId);
                contentStringBuilder.AppendLine(String.Format("{0}  {1}  {2}", this.SafetyFactorTimeserie.Entries.Count(), 2, 1));
                for (int safetyFactorIndex = 0;
                     safetyFactorIndex < this.SafetyFactorTimeserie.Entries.Count();
                     safetyFactorIndex++)
                {
                    contentStringBuilder.AppendLine(String.Format("{0:e}  {1:e}",
                        (this.SafetyFactorTimeserie.Entries[safetyFactorIndex].DateTime - this.SafetyFactorTimeserie.StartDateTime).TotalDays,
                        this.SafetyFactorTimeserie.Entries[safetyFactorIndex].Value));
                }

                string contentString = contentStringBuilder.ToString();
                System.IO.File.WriteAllText(filenameBase + DupuitCalculatorFileIO.safetyFactorFilenameExt, contentString);
            }
        }

        /// <summary>
        /// Write precipitation to file
        /// </summary>
        /// <param name="filenameBase"></param>
        private void WritePrecipitation(string filenameBase)
        {
            var contentStringBuilder = new StringBuilder();
            contentStringBuilder.AppendLine(precipitationId);
            contentStringBuilder.AppendLine(" 2    2    1");
            contentStringBuilder.AppendLine("+0.0000e+00 +0.0000e-02");
            contentStringBuilder.AppendLine("+1.0000e-01 +0.0000e-02");

            // TODO: write real content of precipitation file

            string contentString = contentStringBuilder.ToString();
            System.IO.File.WriteAllText(filenameBase + DupuitCalculatorFileIO.precipitationFilenameExt, contentString);
        }

        /// <summary>
        /// Write monitoring to file
        /// </summary>
        /// <param name="filenameBase"></param>
        private void WriteMonitor(string filenameBase)
        {
            var contentStringBuilder = new StringBuilder();
            contentStringBuilder.AppendLine(monitorId);
            contentStringBuilder.AppendLine("nline    0");

            // TODO: write content of monitor file

            string contentString = contentStringBuilder.ToString();
            System.IO.File.WriteAllText(filenameBase + DupuitCalculatorFileIO.monitorFilenameExt, contentString);
        }

        /// <summary>
        /// Return a geometry2dData with reverserd layer order 
        /// </summary>
        /// <param name="geometry2DData"></param>
        /// <returns></returns>
        private Geometry2DData ReverseGeometry2DData(Geometry2DData geometry2DData)
        {
            if (isReverseLayerOrder)
            {
                var newGeometryData = new Geometry2DData();
                for (int layerIndex = geometry2DData.LayerCount - 1; layerIndex >= 0; layerIndex--)
                {
                    newGeometryData.AddLayer(geometry2DData.GetLayer(layerIndex));
                }
                return newGeometryData;
            }
            else
            {
                return geometry2DData;
            }

        }

        /// <summary>
        /// Return a geometry2dData with reverserd layer order 
        /// </summary>
        /// <param name="geometry2DData"></param>
        /// <returns></returns>
        private DupuitPLLinesTimeSerie ReverseDupuitPLLinesTimeSerie(DupuitPLLinesTimeSerie dupuitPLLinesTimeSerie)
        {
            if (isReverseLayerOrder)
            {
                var newDupuitPLLinesTimeSerie = new DupuitPLLinesTimeSerie();
                foreach (DupuitTimeSerieEntry entry in dupuitPLLinesTimeSerie.Entries)
                {
                    var newEntry = new DupuitTimeSerieEntry() { DupuitPlLines = new DupuitPLLines() };
                    newEntry.TimeInSeconds = entry.TimeInSeconds;
                    newEntry.DupuitPlLines.PLLines = new List<PLLine>();
                    for (int pllineIndex = entry.DupuitPlLines.PLLines.Count - 1; pllineIndex >= 0; pllineIndex--)
                    {
                        newEntry.DupuitPlLines.PLLines.Add(entry.DupuitPlLines.PLLines[pllineIndex]);
                    }
                    newDupuitPLLinesTimeSerie.Entries.Add(newEntry);
                }
                return newDupuitPLLinesTimeSerie;
            }
            else
            {
                return dupuitPLLinesTimeSerie;
            }

        }

        private Soil GetDefaultSoil()
        {
            var soil = new Soil()
            {
                BelowPhreaticLevel = 0.0,
                AbovePhreaticLevel = 0.0,
                FrictionAngle = 0.0,
                Cohesion = 0.0,
                PermeabKx = 0.0
            };
            return soil;
        }

        /// <summary>
        /// Write materials to file
        /// </summary>
        /// <param name="filenameBase"></param>
        private void WriteMaterials(string filenameBase)
        {
            if ((adaptedGeometry2DData != null) && (soilList != null))
            {
                var contentStringBuilder = new StringBuilder();
                contentStringBuilder.AppendLine(materialId);
                contentStringBuilder.AppendLine(String.Format("nLayer   {0}", adaptedGeometry2DData.LayerCount));
                contentStringBuilder.AppendLine(
                    String.Format("{0,8}  {1,14}  {2,14}  {3,14}  {4,14}  {5,14}  {6,14}  {7,14}  {8,14}  {9,14}",
                                  "layer", "bottom(m)", "gW(kN/m3)", "gD(kN/m3)", "alpha(m2/kN)", "phi(o)", "c(kN/m2)", "n(-)",
                                  "Kx(m/d)", "Ky(m/d)"));
                for (int layerIndex = 0; layerIndex < adaptedGeometry2DData.LayerCount; layerIndex++)
                {
                    // In the DikeFlow module the top boundary does not have a soil assigned
                    // In DGeoStability the bottom boundary does not have a soil assigned
                    // So in the conversion between DGeoStability and DikeFlow we have 
                    // to accommodate this. The toplayer does not have a soil assigned 
                    // and each layer should use the soil of the layer before
                    Soil soil;
                    if (layerIndex == 0)
                    {
                        soil = GetDefaultSoil();
                    }
                    else
                    {
                        Geometry2DLayer layerForSoil = adaptedGeometry2DData.GetLayer(layerIndex - 1);
                        int soilIndex = soilList.GetSoilIndexByName(layerForSoil.soilName);
                        soil = (soilIndex >= 0) ? soilList.Soils[soilIndex] : GetDefaultSoil();
                    }

                    Geometry2DLayer layer = adaptedGeometry2DData.GetLayer(layerIndex);
                    contentStringBuilder.AppendLine(
                        String.Format(
                            "{0,8}  {1,14:e}  {2,14:e}  {3,14:e}  {4,14:e}  {5,14:e}  {6,14:e}  {7,14:e}  {8,14:e}  {9,14:e}  {10}",
                            layerIndex,
                            layer.boundaryLine.MinZ(),
                            soil.BelowPhreaticLevel,
                            soil.AbovePhreaticLevel,
                            0.0, //"a", 
                            soil.FrictionAngle,
                            soil.Cohesion,
                            0.3, // "n",
                            soil.PermeabKx * cFactor_mPerSec_2_mPerDay,
                            soil.PermeabKx * cFactor_mPerSec_2_mPerDay,
                            LayerType2String(layer.layerType)));
                }
                contentStringBuilder.AppendLine(String.Format("cDitch(d)      {0:e}", +0.0000e+00));
                contentStringBuilder.AppendLine(String.Format("cRiverBed(d)   {0:e}", +0.0000e+00));
                contentStringBuilder.AppendLine(String.Format("cPolder(d)     {0:e}", +0.0000e+00));
                contentStringBuilder.AppendLine(String.Format("cDike(d)       {0:e}", +0.0000e+00));
                contentStringBuilder.AppendLine(String.Format("beta(m2/kN)    {0:e}", +1.0000e+02));

                string contentString = contentStringBuilder.ToString();
                System.IO.File.WriteAllText(filenameBase + DupuitCalculatorFileIO.materialFilenameExt, contentString);
            }
        }
        /// <summary>
        /// Write calculation options to file
        /// </summary>
        /// <param name="filenameBase"></param>
        private void WriteCalculationOptions(string filenameBase)
        {
            if (dupuitCalculationOptions != null)
            {
                var contentStringBuilder = new StringBuilder();
                contentStringBuilder.AppendLine(calculateId);
                contentStringBuilder.AppendLine(String.Format("version     {0}", dupuitCalculationOptions.Version));
                contentStringBuilder.AppendLine(String.Format("nx          {0}", dupuitCalculationOptions.GridCountX));
                contentStringBuilder.AppendLine(String.Format("ny          {0}", dupuitCalculationOptions.GridCountY));
                contentStringBuilder.AppendLine(String.Format("nr          {0}", dupuitCalculationOptions.GridCircleCount));
                contentStringBuilder.AppendLine(String.Format("cutting     {0}", dupuitCalculationOptions.Cutting));
                contentStringBuilder.AppendLine(String.Format("deformation {0}", dupuitCalculationOptions.Deformation));
                contentStringBuilder.AppendLine(String.Format("generation  {0}", dupuitCalculationOptions.Generation));
                contentStringBuilder.AppendLine(String.Format("x0          {0:e}", dupuitCalculationOptions.GridX0));
                contentStringBuilder.AppendLine(String.Format("x1          {0:e}", dupuitCalculationOptions.GridX1));
                contentStringBuilder.AppendLine(String.Format("y0          {0:e}", dupuitCalculationOptions.GridY0));
                contentStringBuilder.AppendLine(String.Format("y1          {0:e}", dupuitCalculationOptions.GridY1));
                contentStringBuilder.AppendLine(String.Format("r0          {0:e}", dupuitCalculationOptions.GridTangentLineLevel1));
                contentStringBuilder.AppendLine(String.Format("r1          {0:e}", dupuitCalculationOptions.GridTangentLineLevel0));
                contentStringBuilder.AppendLine(String.Format("nxdupuit    {0}", dupuitCalculationOptions.DupuitCountX));
                contentStringBuilder.AppendLine(String.Format("ntdupuit    {0}", dupuitCalculationOptions.DupuitCountT));
                contentStringBuilder.AppendLine(String.Format("stabmethod  {0}", (int)dupuitCalculationOptions.StabMethod));
                contentStringBuilder.AppendLine(String.Format("ntstab      {0}", dupuitCalculationOptions.StabilityCountT));
                contentStringBuilder.AppendLine(String.Format("nxstab      {0}", dupuitCalculationOptions.StabilityCountX));
                contentStringBuilder.AppendLine(String.Format("sTime       {0:e}", dupuitCalculationOptions.TimeStart));
                contentStringBuilder.AppendLine(String.Format("eTime       {0:e}", dupuitCalculationOptions.TimeEnd));
                contentStringBuilder.AppendLine(String.Format("sPolder     {0:e}", dupuitCalculationOptions.PolderS));

                string contentString = contentStringBuilder.ToString();
                System.IO.File.WriteAllText(filenameBase + DupuitCalculatorFileIO.calculateFilenameExt, contentString);
            }
        }

        /// <summary>
        /// Add section string to content
        /// Sections are all mandatory, so when a characteristic point is not available 
        /// the last available point will be used as 
        /// </summary>
        /// <param name="contentStringBuilder"></param>
        /// <param name="characteristicPointType"></param>
        /// <param name="idString"></param>
        /// <returns></returns>
        private void AddSectionStringToContent(StringBuilder contentStringBuilder, CharacteristicPointType characteristicPointType, String idString, ref GeometryPoint lastExistingPoint)
        {
            if (surfaceLine.CharacteristicPoints.GetGeometryPoint(characteristicPointType) != null)
            {
                lastExistingPoint = surfaceLine.CharacteristicPoints.GetGeometryPoint(characteristicPointType);
            }
            String sectionString = String.Format("{0, 14}  {1:e}", idString, lastExistingPoint.X);
            contentStringBuilder.AppendLine(sectionString);
        }

        /// <summary>
        /// Write sections to file
        /// </summary>
        /// <param name="filenameBase"></param>
        private void WriteSections(string filenameBase)
        {
            if (surfaceLine != null)
            {
                var contentStringBuilder = new StringBuilder();
                contentStringBuilder.AppendLine(sectionsId);
                GeometryPoint lastExistingPoint = surfaceLine.CharacteristicPoints.GetGeometryPoint(CharacteristicPointType.SurfaceLevelOutside);
                AddSectionStringToContent(contentStringBuilder, CharacteristicPointType.SurfaceLevelOutside, "xLeftBound", ref lastExistingPoint);
                AddSectionStringToContent(contentStringBuilder, CharacteristicPointType.DikeToeAtRiver, "xOuterToe", ref lastExistingPoint);
                AddSectionStringToContent(contentStringBuilder, CharacteristicPointType.DikeTopAtRiver, "xOuterCrest", ref lastExistingPoint);
                AddSectionStringToContent(contentStringBuilder, CharacteristicPointType.DikeTopAtPolder, "xInnerCrest", ref lastExistingPoint);
                AddSectionStringToContent(contentStringBuilder, CharacteristicPointType.ShoulderBaseInside, "xStartBerm", ref lastExistingPoint);
                AddSectionStringToContent(contentStringBuilder, CharacteristicPointType.ShoulderTopInside, "xEndBerm", ref lastExistingPoint);
                AddSectionStringToContent(contentStringBuilder, CharacteristicPointType.DikeToeAtPolder, "xInnerToe", ref lastExistingPoint);
                AddSectionStringToContent(contentStringBuilder, CharacteristicPointType.DitchDikeSide, "xStartDitch", ref lastExistingPoint);
                AddSectionStringToContent(contentStringBuilder, CharacteristicPointType.DitchPolderSide, "xEndDitch", ref lastExistingPoint);
                AddSectionStringToContent(contentStringBuilder, CharacteristicPointType.SurfaceLevelInside, "xRightBound", ref lastExistingPoint);
                contentStringBuilder.AppendLine(String.Format("nMonitor        {0}", 0));
                //contentStringBuilder.AppendLine(String.Format("xMonitor1       {0:e}", +6.000e+01));
                //contentStringBuilder.AppendLine(String.Format("yMonitor1       {0:e}", +1.000e+00));

                string contentString = contentStringBuilder.ToString();
                System.IO.File.WriteAllText(filenameBase + DupuitCalculatorFileIO.sectionsFilenameExt, contentString);
            }
        }

        /// <summary>
        /// Write PLLines to files
        /// </summary>
        /// <param name="filenameBase"></param>
        private void WritePlLines(string filenameBase)
        {
            if (DupuitPlLinesTimeSerie != null)
            {
                for (int timestepIndex = 0; timestepIndex < DupuitPlLinesTimeSerie.Entries.Count; timestepIndex++)
                {
                    var contentStringBuilder = new StringBuilder();
                    contentStringBuilder.AppendLine(piezolinesId);
                    string extension = String.Format("{0}{1}", DupuitCalculatorFileIO.piezoLinesFilesnameExt, timestepIndex.ToString("D3"));
                    DupuitPLLines dupuitPlLines = DupuitPlLinesTimeSerie.Entries[timestepIndex].DupuitPlLines;

                    contentStringBuilder.AppendLine(String.Format(" time  {0:e}  nline   {1}",
                        DupuitPlLinesTimeSerie.Entries[timestepIndex].TimeInSeconds,
                        dupuitPlLines.PLLines.Count));
                    for (int lineIndex = 0; lineIndex < dupuitPlLines.PLLines.Count; lineIndex++)
                    {
                        contentStringBuilder.AppendLine(String.Format(" {0}  {1}   {2}", dupuitPlLines.PLLines[lineIndex].Points.Count, 2, lineIndex + 1));
                        for (int pointIndex = 0; pointIndex < dupuitPlLines.PLLines[lineIndex].Points.Count; pointIndex++)
                        {
                            contentStringBuilder.AppendLine(String.Format("{0:e}   {1:e}",
                                dupuitPlLines.PLLines[lineIndex].Points[pointIndex].X,
                                dupuitPlLines.PLLines[lineIndex].Points[pointIndex].Z));
                        }
                    }
                    string contentString = contentStringBuilder.ToString();
                    System.IO.File.WriteAllText(Path.ChangeExtension(filenameBase, extension), contentString);
                }
            }
        }

        /// <summary>
        /// Convert layertype to string
        /// </summary>
        /// <param name="layerType"></param>
        /// <returns></returns>
        string LayerType2String(LayerType layerType)
        {
            switch (layerType)
            {
                case LayerType.Top: return layerTypeTop;
                case LayerType.Dike: return layerTypeDike;
                case LayerType.Cover: return layerTypeCover;
                case LayerType.Aquifer1: return layerTypeAquifer1;
                case LayerType.Aquifer2: return layerTypeAquifer2;
                case LayerType.Aquitard: return layerTypeAquitard;

            }
            return null;
        }

        /// <summary>
        /// Read safetyfactors
        /// </summary>
        /// <param name="filenameBase"></param>
        private void ReadSafetyFactors(string filenameBase)
        {
            string filename = filenameBase + DupuitCalculatorFileIO.safetyFactorFilenameExt;
            bool fileFound = File.Exists(filename);
            if (fileFound)
            {
                this.SafetyFactorTimeserie = new TimeSerie();
                int lineIndex = 0;
                string content = System.IO.File.ReadAllText(filename);
                string[] lines = Regex.Split(content, Environment.NewLine);
                ThrowIfFalse(lines[lineIndex].Trim().Equals(safetyFactorId), String.Format("Expected file id '{0}', read file id '{1}'", safetyFactorId, lines[lineIndex]));
                lineIndex++;
                string[] values = Regex.Split(lines[lineIndex].Trim(), "[ ]+"); ;
                int safetyFactorCount = int.Parse(values[0], CultureInfo.InvariantCulture);
                for (int safetyFactorIndex = 0; safetyFactorIndex < safetyFactorCount; safetyFactorIndex++)
                {
                    lineIndex++;
                    values = Regex.Split(lines[lineIndex].Trim(), "[ ]+"); ;
                    var timeSerieEntry = new TimeSerieEntry();
                    timeSerieEntry.DateTime.AddDays(Double.Parse(values[0], CultureInfo.InvariantCulture));
                    timeSerieEntry.Value = Double.Parse(values[1], CultureInfo.InvariantCulture);
                    this.SafetyFactorTimeserie.Entries.Add(timeSerieEntry);
                }
            }
        }

        /// <summary>
        /// Read PL-lines for each time step
        /// </summary>
        /// <param name="filenameBase"></param>
        private void ReadPlLines(string filenameBase)
        {
            this.DupuitPlLinesTimeSerie = new DupuitPLLinesTimeSerie();
            bool fileFound;
            int timestepIndex = 0;
            do
            {
                string extension = String.Format("{0}{1}", DupuitCalculatorFileIO.piezoLinesFilesnameExt, timestepIndex.ToString("D3"));

                string filename = filenameBase + extension;
                fileFound = File.Exists(filename);
                if (fileFound)
                {
                    var dupuitTimeSerieEntry = new DupuitTimeSerieEntry();
                    int lineIndex = 0;
                    string content = System.IO.File.ReadAllText(filename);
                    string[] lines = Regex.Split(content, Environment.NewLine);
                    //lines[0].Trim().Equals();
                    ThrowIfFalse(lines[lineIndex].Trim().Equals(piezolinesId), String.Format("Expected file id '{0}', read file id '{1}'", piezolinesId, lines[lineIndex]));
                    lineIndex++;
                    string[] values = Regex.Split(lines[lineIndex].Trim(), "[ ]+"); ;
                    int plLineCount = int.Parse(values[3], CultureInfo.InvariantCulture);
                    dupuitTimeSerieEntry.TimeInSeconds = Double.Parse(values[1], CultureInfo.InvariantCulture);
                    DupuitPLLines dupuitPlLines = new DupuitPLLines();
                    for (int plLineIndex = 0; plLineIndex < plLineCount; plLineIndex++)
                    {
                        var plLine = new PLLine();
                        lineIndex++;
                        values = Regex.Split(lines[lineIndex].Trim(), "[ ]+"); ;
                        int pointCount = int.Parse(values[0], CultureInfo.InvariantCulture);
                        int lineId = int.Parse(values[2], CultureInfo.InvariantCulture);
                        for (int pointIndex = 0; pointIndex < pointCount; pointIndex++)
                        {
                            lineIndex++;
                            var point = new PLLinePoint();
                            values = Regex.Split(lines[lineIndex].Trim(), "[ ]+");
                            point.X = Double.Parse(values[0], CultureInfo.InvariantCulture);
                            point.Z = Double.Parse(values[1], CultureInfo.InvariantCulture);
                            plLine.Points.Add(point);
                        }
                        dupuitPlLines.PLLines.Add(plLine);
                    }
                    dupuitTimeSerieEntry.DupuitPlLines = dupuitPlLines;
                    timestepIndex++;
                    DupuitPlLinesTimeSerie.Entries.Add(dupuitTimeSerieEntry);
                }
            } while (fileFound);

        }
    }
}