diff --git a/.gitignore b/.gitignore
index 2690ac67..db70e840 100644
--- a/.gitignore
+++ b/.gitignore
@@ -1,3 +1,8 @@
+# pyflowline
+data/**/output/
+examples/**/output/
+examples/**/*.png
+
 # Byte-compiled / optimized / DLL files
 __pycache__/
 *.py[cod]
diff --git a/docs/source/application/application.rst b/docs/source/application/application.rst
index 0d1a2422..925744e3 100644
--- a/docs/source/application/application.rst
+++ b/docs/source/application/application.rst
@@ -24,9 +24,9 @@ The example `run_simulation_mpas.py` script import a few packages and functions.
     import os, sys
     from pathlib import Path
     from os.path import realpath
-    from pyflowline.pyflowline_read_model_configuration_file import pyflowline_read_model_configuration_file
+    from pyflowline.configuration.read_configuration_file import pyflowline_read_configuration_file
 
-The `pyflowline_read_model_configuration_file` function reads in a JSON configuration file and loads all the necessary model parameters. 
+The `pyflowline_read_configuration_file` function reads in a JSON configuration file and loads all the necessary model parameters. 
 
 
 ================
@@ -55,7 +55,7 @@ Check the configuration file:
     if os.path.isfile(sFilename_configuration_in):
         pass
     else:
-        print('This configuration does not exist: ', sFilename_configuration_in )
+        print('The domain configuration file does not exist: ',sFilename_configuration_in )
 
 ================
 Step 4
@@ -70,7 +70,7 @@ Set up case information and read the configuration file.
     Date='20220901'
 
 
-    oPyflowline = pyflowline_read_model_configuration_file(sFilename_configuration_in, \
+    oPyflowline = pyflowline_read_configuration_file(sFilename_configuration_in, \
        iCase_index_in=iCase_index, sDate_in=sDate)
     oPyflowline.aBasin[0].dLatitude_outlet_degree=39.462000
     oPyflowline.aBasin[0].dLongitude_outlet_degree=-76.009300
diff --git a/docs/source/data/data.rst b/docs/source/data/data.rst
index 9fd3b2a7..cb16f843 100644
--- a/docs/source/data/data.rst
+++ b/docs/source/data/data.rst
@@ -54,7 +54,7 @@ PyFlowline uses two JSON-format configuration files to manage all input informat
 
 These files serve as the entry point for setting up and running a PyFlowline case. They can exist wherever the user prefers, but PyFlowline uses the paths specified in these files to locate model inputs and write outputs. Model inputs, outputs, and a recommended directory structure are described in the following two sections.
 
-To create a new PyFlowline case, pass the full path to the parent configuration file to the `pyflowline_read_model_configuration_file` function. This will return a PyFlowline object configured with the values from the file, and it can be used to run the model. See the example notebooks for a demonstration.
+To create a new PyFlowline case, pass the full path to the parent configuration file to the `pyflowline_read_configuration_file` function. This will return a PyFlowline object configured with the values from the file, and it can be used to run the model. See the example notebooks for a demonstration.
 
 Note that the "parent" configuration file contains one block of parameter-value pairs that apply to the entire domain. In contrast, the "child" configuration file contains one block of parameter-value pairs for each watershed. A domain with a single watershed will have a single block in the "child" configuration file, while a domain with multiple watersheds will have multiple blocks.
 
diff --git a/examples/susquehanna/run_simulation_hexagon.py b/examples/susquehanna/run_simulation_hexagon.py
index c76f3ebb..2f72fb7a 100644
--- a/examples/susquehanna/run_simulation_hexagon.py
+++ b/examples/susquehanna/run_simulation_hexagon.py
@@ -2,42 +2,33 @@
 from pathlib import Path
 from os.path import realpath
 
-
-
-#===================================
-#set up workspace path
-#===================================
+#%% set up workspace path
 sPath_parent = str(Path(__file__).parents[2]) # data is located two dir's up
 import sys
 sys.path.append(sPath_parent)
-from pyflowline.configuration.pyflowline_read_configuration_file import pyflowline_read_model_configuration_file
+
+from pyflowline.configuration.read_configuration_file import pyflowline_read_configuration_file
+
 
 sPath_data = realpath( sPath_parent +  '/data/susquehanna' )
 sWorkspace_input =  str(Path(sPath_data)  /  'input')
 sWorkspace_output=  str(Path(sPath_data)  /  'output')
 
-#===================================
-#you need to update this file based on your own case study
-#===================================
+#%% you need to update this file based on your own case study
 sFilename_configuration_in = realpath( sPath_parent +  '/data/susquehanna/input/pyflowline_susquehanna_hexagon.json' )
 if os.path.isfile(sFilename_configuration_in):
     pass
 else:
     print('This configuration does not exist: ', sFilename_configuration_in )
 
-#===================================
-#setup case information
-#===================================
+#%% setup case information
 iCase_index = 1
 dResolution_meter = 50000
 
 sMesh = 'hexagon'
 sDate='20230101'
 
-
-
-
-oPyflowline = pyflowline_read_model_configuration_file(sFilename_configuration_in, \
+oPyflowline = pyflowline_read_configuration_file(sFilename_configuration_in, \
 iCase_index_in=iCase_index, dResolution_meter_in=dResolution_meter, sDate_in=sDate)
 oPyflowline.aBasin[0].dLatitude_outlet_degree=39.462000
 oPyflowline.aBasin[0].dLongitude_outlet_degree=-76.009300
@@ -50,4 +41,4 @@
 
 print('Finished')
 
-
+# %%
diff --git a/examples/susquehanna/run_simulation_latlon.py b/examples/susquehanna/run_simulation_latlon.py
index 590b566d..db1013bb 100644
--- a/examples/susquehanna/run_simulation_latlon.py
+++ b/examples/susquehanna/run_simulation_latlon.py
@@ -2,40 +2,36 @@
 from pathlib import Path
 from os.path import realpath
 
-#===================================
-#set up workspace path
-#===================================
+#%% set up workspace path
 sPath_parent = str(Path(__file__).parents[2]) # data is located two dir's up
 sPath_parent = str(Path(__file__).parents[2]) # data is located two dir's up
 import sys
 sys.path.append(sPath_parent)
-from pyflowline.configuration.pyflowline_read_configuration_file import pyflowline_read_model_configuration_file
+
+from pyflowline.configuration.pyflowline_read_configuration_file import pyflowline_read_configuration_file
+
 sPath_data = realpath( sPath_parent +  '/data/susquehanna' )
 sWorkspace_input =  str(Path(sPath_data)  /  'input')
 sWorkspace_output=  str(Path(sPath_data)  /  'output')
 
-#===================================
-#you need to update this file based on your own case study
-#===================================
-sFilename_configuration_in = realpath( sPath_parent +  '//data/susquehanna/input//pyflowline_susquehanna_latlon.json' )
+#%% you need to update this file based on your own case study
+sFilename_configuration_in = realpath( sPath_parent + '//data/susquehanna/input//pyflowline_susquehanna_latlon.json' )
 if os.path.isfile(sFilename_configuration_in):
     pass
 else:
     print('This configuration does not exist: ', sFilename_configuration_in )
 
-#===================================
-#setup case information
-#===================================
+#%% setup case information
 iCase_index = 3
 dResolution_meter = 50000
 
 sMesh = 'latlon'
 sDate='20230101'
 
-
-
-oPyflowline = pyflowline_read_model_configuration_file(sFilename_configuration_in, \
+oPyflowline = pyflowline_read_configuration_file(sFilename_configuration_in, \
 iCase_index_in=iCase_index, dResolution_meter_in=dResolution_meter, sDate_in=sDate)
+
+#%%
 oPyflowline.aBasin[0].dLatitude_outlet_degree=39.462000
 oPyflowline.aBasin[0].dLongitude_outlet_degree=-76.009300
 oPyflowline.setup()
diff --git a/examples/susquehanna/run_simulation_mpas.py b/examples/susquehanna/run_simulation_mpas.py
index 00ce5396..15315023 100644
--- a/examples/susquehanna/run_simulation_mpas.py
+++ b/examples/susquehanna/run_simulation_mpas.py
@@ -1,98 +1,175 @@
 import os, sys
 from pathlib import Path
-from os.path import realpath
-#===================================
-#set up workspace path
-#===================================
-sPath_parent = str(Path(__file__).parents[2]) # data is located two dir's up
-
-sys.path.append(sPath_parent)
-from pyflowline.configuration.change_json_key_value import change_json_key_value
-from pyflowline.configuration.pyflowline_read_configuration_file import pyflowline_read_model_configuration_file
-
-sPath_data = realpath( sPath_parent +  '/data/susquehanna' )
-sWorkspace_input =  str(Path(sPath_data)  /  'input')
-sWorkspace_output=  str(Path(sPath_data)  /  'output')
-
-
-#===================================
-#you need to update this file based on your own case study
-#===================================
-sFilename_configuration_in = realpath( sPath_parent +  '/data/susquehanna/input/pyflowline_susquehanna_mpas.json' )
-if os.path.isfile(sFilename_configuration_in):
-    pass
-else:
-    print('This configuration does not exist: ', sFilename_configuration_in )
 
-#===================================
-#setup case information
-#===================================
+
+from pyflowline.configuration.read_configuration_file import pyflowline_read_configuration_file
+from pyflowline.configuration.change_json_key_value import pyflowline_change_json_key_value
+from pyflowline.configuration import path_manager as pyflowline_path_manager
+
+
+#%% Define the case information (configuration file parameters)
+sDomainName = 'susquehanna'
 iCase_index = 1
 iFlag_simulation = 0
 iFlag_visualization = 1
 sMesh = 'mpas'
-sDate='20230701'
+sDate = '20230701'
+
+#%% Define workspace paths and configuration file path parameters
+oPath_parent = pyflowline_path_manager.pyflowline_project_root()
+sys.path.append(str(oPath_parent))
+
+# Define the full path to the input and output folders.
+oFolder_input = oPath_parent.joinpath(
+    'data', sDomainName, 'input')
+oFolder_output = oPath_parent.joinpath(
+    'data', sDomainName, 'output')
+
+# Define the full path to the domain ("parent") configuration file.
+oFilename_domain_config = oFolder_input.joinpath(
+    'pyflowline_susquehanna_mpas.json')
 
-sFolder_data = os.path.join(sPath_parent, 'data')
-sFolder_data_susquehanna =  os.path.join(sFolder_data, 'susquehanna')
-sFolder_input = os.path.join(sFolder_data_susquehanna, 'input')
-sFilename_flowline = realpath( os.path.join(sFolder_input, 'flowline.geojson') )
-sFilename_basins =  realpath( os.path.join(sFolder_input , 'pyflowline_susquehanna_basins.json' ))
-change_json_key_value(sFilename_basins, 'sFilename_flowline_filter', sFilename_flowline, iFlag_basin_in=1)
+# Define the full path to the individual basin ("child") configuration file.
+oFilename_basins_config = oFolder_input.joinpath(
+    'pyflowline_susquehanna_basins.json')
 
-oPyflowline = pyflowline_read_model_configuration_file(sFilename_configuration_in,
-    iCase_index_in=iCase_index, sDate_in=sDate)
+# Define the full path to the MPAS mesh file
+oFilename_mesh_netcdf = oFolder_input.joinpath(
+    'lnd_cull_mesh.nc')
 
-#take a look at the model parameters
+# Define the full path to the input flowline file.
+oFilename_flowline = oFolder_input.joinpath(
+    'flowline.geojson')
 
+# Define the full path to the boundary file used to clip the mesh.
+oFilename_mesh_boundary = oFolder_input.joinpath(
+    'boundary_wgs.geojson')
+
+# Confirm that the domain configuration file exists.
+if os.path.isfile(oFilename_domain_config):
+    pass
+else:
+    print('The domain configuration file does not exist: ', oFilename_domain_config)
+
+#%% Update the domain (parent) configuration file
+
+# The following parameters need to be set in the configuration file:
+#   sWorkspace_output: full/path/to/output
+#   sFilename_mesh_netcdf: full/path/to/lnd_cull_mesh.nc
+#   sFilename_mesh_boundary: full/path/to/boundary_wgs.geojson
+# 	sFilename_basins: full/path/to/pyflowline_susquehanna_basins.json.
+
+# The json file will be overwritten, you may want to make a copy of it first.
+
+# Set the path to the output folder
+# Pass the configuration filename followed by a single key-value pair.
+pyflowline_change_json_key_value(
+    oFilename_domain_config,
+    'sWorkspace_output', str(oFolder_output))
+
+# Set the path to the mpas mesh file
+pyflowline_change_json_key_value(
+    oFilename_domain_config,
+    'sFilename_mesh_netcdf', str(oFilename_mesh_netcdf))
+
+# Set the path to the boundary file used to clip the mesh
+pyflowline_change_json_key_value(
+    oFilename_domain_config,
+    'sFilename_mesh_boundary', str(oFilename_mesh_boundary))
+
+# Set the path to the basin ("child") configuration file
+pyflowline_change_json_key_value(
+    oFilename_domain_config,
+    'sFilename_basins', str(oFilename_basins_config))
+
+# These parameters are not strictly necessary, but will reduce the potential for errors or confusion.
+pyflowline_change_json_key_value(
+    oFilename_domain_config,
+    'sFilename_model_configuration', str(oFilename_domain_config))
+
+pyflowline_change_json_key_value(
+    oFilename_domain_config,
+    'sWorkspace_data', str(oPath_parent.joinpath('data')))
+
+#%% Update the basin ("child") configuration file
+
+# Set the path to the flowline file.
+pyflowline_change_json_key_value(
+    oFilename_basins_config,
+    'sFilename_flowline_filter', str(oFilename_flowline),
+    iFlag_basin_in=1) # Set iFlag_basin_in=1 when changing the basin configuration file.
+
+#%% Read the configuration file
+oPyflowline = pyflowline_read_configuration_file(
+    oFilename_domain_config,
+    iCase_index_in=iCase_index,
+    sDate_in=sDate)
+
+#%% Check the model parameters
 oPyflowline.pyflowline_print()
 
-#now we can change the following model parameters
-#there are two ways to change the model parameters
-#use a function or assign a value directly
+#%% Now we can change some model parameters
 
-oPyflowline.pyflowline_change_model_parameter(sVariable_in='sWorkspace_output', sValue_in=sWorkspace_output)
+# There are two ways to change the model parameters: 1) use a function, or 2) assign a value directly. Use the function:
+oPyflowline.pyflowline_change_model_parameter(
+    'sWorkspace_output', str(oFolder_output))
 
-#if you need to change a parameter for a basin instead of the whole model domain, use the iFlag_basin_in option, this will change all the basins
-oPyflowline.pyflowline_change_model_parameter(sVariable_in='dLatitude_outlet_degree', sValue_in=39.462000, iFlag_basin_in=1)
-oPyflowline.pyflowline_change_model_parameter(sVariable_in='dLongitude_outlet_degree', sValue_in=-76.009300, iFlag_basin_in=1)
+# To change a parameter for a basin instead of the whole model domain, use the iFlag_basin_in option, this will change the value for all of the basins in the basin configuration file.
+oPyflowline.pyflowline_change_model_parameter(
+    'dLatitude_outlet_degree', 39.462000,
+    iFlag_basin_in=1)
 
+oPyflowline.pyflowline_change_model_parameter(
+    'dLongitude_outlet_degree', -76.009300,
+    iFlag_basin_in=1)
 
-#the second way is to assign a value directly
+# The second way is to assign a value directly
 oPyflowline.aBasin[0].dLatitude_outlet_degree=39.462000
 oPyflowline.aBasin[0].dLongitude_outlet_degree=-76.009300
 
-#oPyflowline.setup()
-if iFlag_visualization ==1:
-    #oPyflowline.plot(sVariable_in = 'flowline_filter', sFilename_output_in = 'filter_flowline.png'  )
+#%%
+oPyflowline.pyflowline_setup()
+
+#%%
+if iFlag_visualization == 1:
+    oPyflowline.plot(sVariable_in='flowline_filter', sFilename_output_in='filter_flowline.png')
     pass
 
+#%%
 if iFlag_simulation == 1:
-    #oPyflowline.pyflowline_flowline_simplification()
+    oPyflowline.pyflowline_flowline_simplification()
     pass
 
+#%%
 if iFlag_visualization == 1:
 
-    aExtent_meander = [-76.5,-76.2, 41.6,41.9]
-    #oPyflowline.plot( sVariable_in='flowline_simplified' , sFilename_output_in = 'flowline_simplified.png' )
-    #oPyflowline.plot( sVariable_in='flowline_simplified' , sFilename_output_in = 'flowline_simplified_zoom.png', aExtent_in =aExtent_meander )
-
+    aExtent_meander = [-76.5, -76.2, 41.6, 41.9]
+    # oPyflowline.plot(sVariable_in='flowline_simplified', sFilename_output_in='flowline_simplified.png')
+    # oPyflowline.plot(sVariable_in='flowline_simplified',
+    #                  sFilename_output_in='flowline_simplified_zoom.png',
+    #                  aExtent_in=aExtent_meander)
     pass
 
+#%%
 if iFlag_simulation == 1:
     aCell = oPyflowline.pyflowline_mesh_generation()
 
+#%%
 if iFlag_visualization == 1:
-    oPyflowline.plot( sVariable_in='mesh', sFilename_output_in = 'mesh.png' )
+    # This will fail unless iFlag_simulation is true or has been run
+    # oPyflowline.plot(sVariable_in='mesh', sFilename_output_in='mesh.png' )
     pass
-
+#%%
 if iFlag_simulation == 1:
     oPyflowline.pyflowline_reconstruct_topological_relationship(aCell)
 
+#%%
 if iFlag_visualization == 1:
-    oPyflowline.plot(  sVariable_in='overlap', sFilename_output_in = 'mesh_w_flowline.png',)
+    # oPyflowline.plot(sVariable_in='overlap', sFilename_output_in='mesh_w_flowline.png')
     pass
 
-oPyflowline.export()
+# oPyflowline.pyflowline_export()
 
 print('Finished')
+
+# %%
diff --git a/examples/susquehanna/run_simulation_square.py b/examples/susquehanna/run_simulation_square.py
index 0666d792..2feeb24e 100644
--- a/examples/susquehanna/run_simulation_square.py
+++ b/examples/susquehanna/run_simulation_square.py
@@ -2,39 +2,30 @@
 from pathlib import Path
 from os.path import realpath
 
-
-#===================================
-#set up workspace path
-#===================================
+#%% set up workspace path
 sPath_parent = str(Path(__file__).parents[2]) # data is located two dir's up
 import sys
 sys.path.append(sPath_parent)
-from pyflowline.configuration.pyflowline_read_configuration_file import pyflowline_read_model_configuration_file
+from pyflowline.configuration.read_configuration_file import pyflowline_read_configuration_file
 sPath_data = realpath( sPath_parent +  '/data/susquehanna' )
 sWorkspace_input =  str(Path(sPath_data)  /  'input')
 sWorkspace_output=  str(Path(sPath_data)  /  'output')
 
-#===================================
-#you need to update this file based on your own case study
-#===================================
-sFilename_configuration_in = realpath( sPath_parent +  '/data/susquehanna/input//pyflowline_susquehanna_square.json' )
+#%% you need to update this file based on your own case study
+sFilename_configuration_in = realpath( sPath_parent +  '/data/susquehanna/input/pyflowline_susquehanna_square.json' )
 if os.path.isfile(sFilename_configuration_in):
     pass
 else:
     print('This configuration does not exist: ', sFilename_configuration_in )
 
-#===================================
-#setup case information
-#===================================
+#%% setup case information
 iCase_index = 2
 dResolution_meter =  50000
 
 sMesh = 'square'
 sDate='20230101'
 
-
-
-oPyflowline = pyflowline_read_model_configuration_file(sFilename_configuration_in, \
+oPyflowline = pyflowline_read_configuration_file(sFilename_configuration_in, \
 iCase_index_in=iCase_index, dResolution_meter_in=dResolution_meter, sDate_in=sDate)
 oPyflowline.aBasin[0].dLatitude_outlet_degree=39.462000
 oPyflowline.aBasin[0].dLongitude_outlet_degree=-76.009300
diff --git a/notebooks/mpas_example.ipynb b/notebooks/mpas_example.ipynb
index 982c8734..80ebc15d 100644
--- a/notebooks/mpas_example.ipynb
+++ b/notebooks/mpas_example.ipynb
@@ -41,80 +41,8 @@
     "from pathlib import Path\n",
     "from os.path import realpath\n",
     "import importlib.util\n",
-    "#check dependencies\n",
-    "\n",
-    "iFlag_numpy = importlib.util.find_spec(\"numpy\")\n",
-    "iFlag_gdal = importlib.util.find_spec(\"osgeo\")\n",
-    "iFlag_netcdf4 = importlib.util.find_spec(\"netCDF4\")\n",
-    "iFlag_geopandas = importlib.util.find_spec(\"geopandas\")\n",
-    "iFlag_requests = importlib.util.find_spec(\"requests\") #this one is only for downloading the data"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": 2,
-   "id": "4b2ca848",
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "if iFlag_numpy is not None:\n",
-    "    pass\n",
-    "else:\n",
-    "    print(\"numpy is not installed\")"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": 3,
-   "id": "b4ba011d",
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "if iFlag_gdal is not None:\n",
-    "    pass\n",
-    "else:\n",
-    "    print(\"gdal is not installed\")"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": 4,
-   "id": "f5c78ad3",
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "if iFlag_netcdf4 is not None:\n",
-    "    pass\n",
-    "else:\n",
-    "    print(\"netCDF4 is not installed\")"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": 5,
-   "id": "4d2b0424",
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "if iFlag_geopandas is not None:\n",
-    "    pass\n",
-    "else:\n",
-    "    print(\"geopandas is not installed\")"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": 6,
-   "id": "21dd10ee",
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "if iFlag_requests is not None:\n",
-    "    pass\n",
-    "else:\n",
-    "    print(\"requests is not installed\")\n",
-    "\n",
-    "import requests"
+    "import requests\n",
+    "#check dependencies\n"
    ]
   },
   {
@@ -143,17 +71,15 @@
    "source": [
     "#now add the pyflowline into the Python path.\n",
     "sPath_notebook = Path().resolve()\n",
-    "sPath_parent = str(Path().resolve().parents[0]) \n",
+    "sPath_parent = str(Path().resolve().parents[0])\n",
     "print(sPath_parent)\n",
     "#check pyflowline installation\n",
-    "iFlag_pyflowline = importlib.util.find_spec(\"pyflowline\") \n",
+    "iFlag_pyflowline = importlib.util.find_spec(\"pyflowline\")\n",
     "if iFlag_pyflowline is not None:\n",
     "    pass\n",
     "else:\n",
     "    print('The pyflowline package is not installed. We will use the current path to set it up.')\n",
-    "    sys.path.append(sPath_parent)\n",
-    "\n",
-    "\n"
+    "    sys.path.append(sPath_parent)"
    ]
   },
   {
@@ -183,7 +109,6 @@
     }
    ],
    "source": [
-    "\n",
     "#download the MPAS mesh from the github release\n",
     "sFilename_mpas = 'https://github.com/changliao1025/pyflowline/releases/download/0.2.0/lnd_cull_mesh.nc'\n",
     "\n",
@@ -193,6 +118,7 @@
     "sFolder_input = os.path.join(sFolder_data_susquehanna, 'input')\n",
     "\n",
     "sFilename_download = os.path.join(sFolder_input, 'mpas_mesh.nc')\n",
+    "\n",
     "#check whether the file exists\n",
     "if os.path.exists(sFilename_download):\n",
     "    print(f\"File '{sFilename_download}' exists. Skip downloading.\")\n",
@@ -229,8 +155,8 @@
    "source": [
     "#step 3\n",
     "#load the read configuration function\n",
-    "from pyflowline.change_json_key_value import change_json_key_value\n",
-    "from pyflowline.pyflowline_read_model_configuration_file import pyflowline_read_model_configuration_file"
+    "from pyflowline.configuration.change_json_key_value import pyflowline_change_json_key_value\n",
+    "from pyflowline.configuration.read_configuration_file import pyflowline_read_configuration_file"
    ]
   },
   {
@@ -307,7 +233,6 @@
     }
    ],
    "source": [
-    "\n",
     "#we can check what is the content of this json file\n",
     "import json\n",
     "with open(sFilename_configuration_in, 'r') as pJSON:\n",
@@ -362,17 +287,38 @@
    "source": [
     "#we need to update a few parameters in the configuration file before we can create the flowline object\n",
     "#the json file will be overwritten, you may want to make a copy of it first\n",
-    "sFilename_basins =  realpath( os.path.join(sFolder_input , 'pyflowline_susquehanna_basins.json' ))\n",
-    "sFilename_mesh_boundary = realpath(os.path.join(sFolder_input, 'boundary_wgs.geojson')) #boundary to clip mesh\n",
+    "sFilename_basins = realpath(os.path.join(\n",
+    "    sFolder_input, 'pyflowline_susquehanna_basins.json'))\n",
+    "\n",
+    "sFilename_mesh_boundary = realpath(os.path.join(\n",
+    "    sFolder_input, 'boundary_wgs.geojson')) #boundary to clip mesh\n",
+    "\n",
     "sWorkspace_output = os.path.join(sFolder_data_susquehanna, 'output')\n",
-    "change_json_key_value(sFilename_configuration_in, 'sFilename_mesh_netcdf', sFilename_download) #mpas file we just downloaded\n",
-    "change_json_key_value(sFilename_configuration_in, 'sFilename_mesh_boundary', sFilename_mesh_boundary) \n",
-    "change_json_key_value(sFilename_configuration_in, 'sFilename_basins', sFilename_basins) #individual basin configuration file\n",
-    "change_json_key_value(sFilename_configuration_in, 'sWorkspace_output', sWorkspace_output) #output folder\n",
+    "\n",
+    "pyflowline_change_json_key_value(sFilename_configuration_in,\n",
+    "                                 'sFilename_mesh_netcdf',\n",
+    "                                 sFilename_download) #mpas file we just downloaded\n",
+    "\n",
+    "pyflowline_change_json_key_value(sFilename_configuration_in,\n",
+    "                                 'sFilename_mesh_boundary',\n",
+    "                                 sFilename_mesh_boundary)\n",
+    "\n",
+    "pyflowline_change_json_key_value(sFilename_configuration_in,\n",
+    "                                 'sFilename_basins',\n",
+    "                                 sFilename_basins) #individual basin configuration file\n",
+    "\n",
+    "pyflowline_change_json_key_value(sFilename_configuration_in,\n",
+    "                                 'sWorkspace_output',\n",
+    "                                 sWorkspace_output) #output folder\n",
     "\n",
     "#now change basin configuration file\n",
-    "sFilename_flowline = realpath( os.path.join(sFolder_input, 'flowline.geojson') )\n",
-    "change_json_key_value(sFilename_basins, 'sFilename_flowline_filter', sFilename_flowline, iFlag_basin_in=1) #user provided flowline\n"
+    "sFilename_flowline = realpath(os.path.join(\n",
+    "    sFolder_input, 'flowline.geojson'))\n",
+    "\n",
+    "pyflowline_change_json_key_value(sFilename_basins,\n",
+    "                                 'sFilename_flowline_filter',\n",
+    "                                 sFilename_flowline,\n",
+    "                                 iFlag_basin_in=1) #user provided flowline"
    ]
   },
   {
@@ -404,8 +350,11 @@
     "#iCase_index_in: this is an ID to identify the simulation case\n",
     "#sMesh_type_in: this specifies the mesh type ('mpas' in this example)\n",
     "#sDate_in: this specifies the date of the simulation, the final output folder will have a pattern such as 'pyflowline20230901001', where pyflowline is model, 20230901 is the date, and 001 is the case index.\n",
-    "oPyflowline = pyflowline_read_model_configuration_file(sFilename_configuration_in, iCase_index_in=iCase_index, \n",
-    "                                                    sMesh_type_in= sMesh_type, sDate_in=sDate)"
+    "oPyflowline = pyflowline_read_configuration_file(\n",
+    "    sFilename_configuration_in,\n",
+    "    iCase_index_in=iCase_index,\n",
+    "    sMesh_type_in=sMesh_type,\n",
+    "    sDate_in=sDate)"
    ]
   },
   {
@@ -419,14 +368,19 @@
     "#only a list of parameters can be changed, for full list, please check the documentation\n",
     "#in this example, we will change the mesh file name\n",
     "#the function will check the data type, if incorrect data type is provided, it will raise an error\n",
-    "oPyflowline.pyflowline_change_model_parameter('sFilename_mesh_netcdf', sFilename_download) #because the mpas mesh already contains elevation, we do not need to set the elevation file name\n",
-    "oPyflowline.pyflowline_change_model_parameter('sFilename_mesh_boundary', sFilename_mesh_boundary)\n",
+    "oPyflowline.pyflowline_change_model_parameter(\n",
+    "    'sFilename_mesh_netcdf', sFilename_download)\n",
+    "#because the mpas mesh already contains elevation, we do not need to set the elevation file name\n",
+    "oPyflowline.pyflowline_change_model_parameter(\n",
+    "    'sFilename_mesh_boundary', sFilename_mesh_boundary)\n",
     "\n",
     "#we can also set for individual basin in the domain, in this example, we only has one basin.\n",
     "#remember that, each basin can have different parameters, so if you want to set them different (for example, basin 1 has no dam, but basin 2 has dam), you should edit the basin json instead using this function.\n",
     "#because change_model_parameter will set all the basin using the same parameter in current version\n",
     "#must set iFlag_basin_in = 1 for basin parameter\n",
-    "oPyflowline.pyflowline_change_model_parameter('iFlag_dam', 0, iFlag_basin_in= 1)\n"
+    "oPyflowline.pyflowline_change_model_parameter(\n",
+    "    'iFlag_dam', 0,\n",
+    "    iFlag_basin_in= 1)"
    ]
   },
   {
@@ -515,8 +469,12 @@
    "source": [
     "#another important setting for basin is the approximate outlet location\n",
     "#you can set it using the change_model_parameter function\n",
-    "oPyflowline.pyflowline_change_model_parameter('dLongitude_outlet_degree', -76.0093, iFlag_basin_in = 1)\n",
-    "oPyflowline.pyflowline_change_model_parameter('dLatitude_outlet_degree', 39.4620, iFlag_basin_in = 1)"
+    "oPyflowline.pyflowline_change_model_parameter(\n",
+    "    'dLongitude_outlet_degree', -76.0093,\n",
+    "    iFlag_basin_in = 1)\n",
+    "oPyflowline.pyflowline_change_model_parameter(\n",
+    "    'dLatitude_outlet_degree', 39.4620,\n",
+    "    iFlag_basin_in = 1)"
    ]
   },
   {
@@ -614,7 +572,7 @@
     }
    ],
    "source": [
-    "#setup the model       \n",
+    "#setup the model\n",
     "oPyflowline.pyflowline_setup()"
    ]
   },
@@ -646,20 +604,14 @@
    ],
    "source": [
     "\n",
-    "\n",
-    "if iFlag_geopandas is not None:\n",
-    "    import geopandas as gpd\n",
-    "    import matplotlib.pyplot as plt\n",
-    "    #use the geopanda package\n",
-    "    #the raw/original geojson file \n",
-    "    sFilename_geojson = oPyflowline.aBasin[0].sFilename_flowline_filter_geojson\n",
-    "    gdf = gpd.read_file(sFilename_geojson)\n",
-    "    gdf.plot()\n",
-    "    plt.show()\n",
-    "else:\n",
-    "    print('The visulization packages are not installed.')\n",
-    "pass\n",
-    "    "
+    "import geopandas as gpd\n",
+    "import matplotlib.pyplot as plt\n",
+    "#use the geopanda package\n",
+    "#the raw/original geojson file\n",
+    "sFilename_geojson = oPyflowline.aBasin[0].sFilename_flowline_filter_geojson\n",
+    "gdf = gpd.read_file(sFilename_geojson)\n",
+    "gdf.plot()\n",
+    "plt.show()\n"
    ]
   },
   {
@@ -762,8 +714,7 @@
     "sFilename_geojson = oPyflowline.aBasin[0].sFilename_flowline_simplified\n",
     "gdf = gpd.read_file(sFilename_geojson)\n",
     "gdf.plot()\n",
-    "plt.show()\n",
-    "pass"
+    "plt.show()\n"
    ]
   },
   {
@@ -1425,13 +1376,11 @@
     }
    ],
    "source": [
-    "\n",
     "#visualize the mesh\n",
     "sFilename_geojson = oPyflowline.sFilename_mesh\n",
     "gdf = gpd.read_file(sFilename_geojson)\n",
     "gdf.plot()\n",
-    "plt.show()\n",
-    "pass"
+    "plt.show()\n"
    ]
   },
   {
@@ -1575,8 +1524,6 @@
     }
    ],
    "source": [
-    "\n",
-    "\n",
     "#plot both the mesh and the flowline\n",
     "file1_path = oPyflowline.sFilename_mesh\n",
     "file2_path = oPyflowline.aBasin[0].sFilename_flowline_conceptual\n",
@@ -1585,8 +1532,7 @@
     "fig, ax = plt.subplots()\n",
     "gdf1.plot(ax=ax, color='blue')\n",
     "gdf2.plot(ax=ax, color='red')\n",
-    "plt.show()\n",
-    "pass"
+    "plt.show()\n"
    ]
   },
   {
@@ -1721,8 +1667,7 @@
    "source": [
     "with open(oPyflowline.sFilename_mesh_info, 'r') as pJSON:\n",
     "    parsed = json.load(pJSON)\n",
-    "    print(json.dumps(parsed[0], indent=4))\n",
-    "    "
+    "    print(json.dumps(parsed[0], indent=4))"
    ]
   },
   {
@@ -1762,8 +1707,10 @@
     }
    ],
    "source": [
-    "sFilename_flowline_conceptual_info= os.path.join(str(Path(oPyflowline.aBasin[0].sWorkspace_output_basin)  ), oPyflowline.aBasin[0].sFilename_flowline_conceptual_info  )    \n",
-    "                \n",
+    "sFilename_flowline_conceptual_info = os.path.join(\n",
+    "    str(Path(oPyflowline.aBasin[0].sWorkspace_output_basin)),\n",
+    "    oPyflowline.aBasin[0].sFilename_flowline_conceptual_info)\n",
+    "\n",
     "with open(sFilename_flowline_conceptual_info, 'r') as pJSON:\n",
     "    parsed = json.load(pJSON)\n",
     "    print(json.dumps(parsed[0], indent=4))"
@@ -1804,7 +1751,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.11.5"
+   "version": "3.9.19"
   },
   "vscode": {
    "interpreter": {
diff --git a/pyflowline/classes/_hpc.py b/pyflowline/classes/_hpc.py
index 79ef6a54..e52fa7a8 100644
--- a/pyflowline/classes/_hpc.py
+++ b/pyflowline/classes/_hpc.py
@@ -16,12 +16,12 @@ def _pyflowline_create_hpc_job(self, sSlurm_in=None):
 
     sLine = '#!/qfs/people/liao313/.conda/envs/pyflowline/bin/' + 'python3' + '\n'
     ofs_pyflowline.write(sLine)
-    sLine = 'from pyflowline.pyflowline_read_model_configuration_file import pyflowline_read_model_configuration_file' + '\n'
+    sLine = 'from pyflowline.configuration.read_configuration_file import pyflowline_read_configuration_file' + '\n'
     ofs_pyflowline.write(sLine)
     sLine = 'sFilename_configuration_in = ' + '"' + \
         self.sFilename_model_configuration + '"\n'
     ofs_pyflowline.write(sLine)
-    sLine = 'oPyflowline = pyflowline_read_model_configuration_file(sFilename_configuration_in,'\
+    sLine = 'oPyflowline = pyflowline_read_configuration_file(sFilename_configuration_in,'\
         + 'iCase_index_in=' + str(self.iCase_index) + ','\
         + 'iResolution_index_in=' + str(self.iResolution_index) + ','\
         + 'dResolution_meter_in=' + "{:0f}".format(self.dResolution_meter) + ','\
diff --git a/pyflowline/configuration/change_json_key_value.py b/pyflowline/configuration/change_json_key_value.py
index 46da235f..13b3737f 100644
--- a/pyflowline/configuration/change_json_key_value.py
+++ b/pyflowline/configuration/change_json_key_value.py
@@ -1,16 +1,39 @@
 import json
+from pathlib import Path
+import shutil
+import tempfile
 
-def change_json_key_value(sFilename_json_in, sKey, new_value, iFlag_basin_in= None):
-    # Read the JSON file
+def pyflowline_change_json_key_value(sFilename_json_in, sKey, new_value, iFlag_basin_in=None):
+    """
+    Change the value associated with the specified key in a JSON file safely.
+
+    Args:
+        sFilename_json_in (str or Path): Path to the JSON file.
+        sKey (str): Key whose value needs to be changed.
+        new_value: New value to set for the specified key.
+        iFlag_basin_in (bool, optional): Flag indicating if the key is for a basin configuration. Defaults to None.
+    """
+    # Ensure input filename is a string
+    sFilename_json_in = str(sFilename_json_in)
+
+    # Convert new_value to string if it's a Path object
+    if isinstance(new_value, Path):
+        new_value = str(new_value)
+
+    # Read the original JSON data
     with open(sFilename_json_in, 'r') as file:
         data = json.load(file)
 
-    # Update the value associated with the specified sKey
+    # Update the value associated with the specified key
     if iFlag_basin_in is None:
         data[sKey] = new_value
     else:
         data[0][sKey] = new_value
 
-    # Write the updated data back to the JSON file
-    with open(sFilename_json_in, 'w') as file:
-        json.dump(data, file, indent=4)
\ No newline at end of file
+    # Write the updated data to a temporary file
+    with tempfile.NamedTemporaryFile(delete=False, mode='w', suffix='.json') as temp_file:
+        json.dump(data, temp_file, indent=4)
+        temp_path = temp_file.name
+
+    # Replace the original file with the updated temporary file
+    shutil.move(temp_path, sFilename_json_in)
diff --git a/pyflowline/configuration/create_template_configuration_file.py b/pyflowline/configuration/create_template_configuration_file.py
new file mode 100644
index 00000000..81ad2bc7
--- /dev/null
+++ b/pyflowline/configuration/create_template_configuration_file.py
@@ -0,0 +1,187 @@
+import os
+from pathlib import Path
+#use this function to generate an initial json file for hexwatershed
+import json
+#once it's generated, you can modify it and use it for different simulations
+from pyflowline.classes.pycase import flowlinecase
+from pyflowline.classes.basin import pybasin
+
+def pyflowline_create_template_basin_configuration_file(
+        sFilename_basins_json,
+        nBasin,
+        sWorkspace_input_in,
+        sWorkspace_output_in):
+    """generate basin configuration
+
+    Args:
+        sFilename_basins_json (str or Path): the filename
+        nBasin (int): the total number of basin
+        sWorkspace_input_in (str or Path): the input data path
+        sWorkspace_output (str or Path): the output path
+
+    Returns:
+        basin: a basin object
+    """
+
+    # Ensure input pathnames are strings
+    sFilename_basins_json = str(sFilename_basins_json)
+    sWorkspace_input_in = str(sWorkspace_input_in)
+    sWorkspace_output_in = str(sWorkspace_output_in)
+
+    aBasin_out = list()
+    for i in range(nBasin):
+        sBasin =  "{:03d}".format(i+1)   
+        aConfig_basin = {}
+        aConfig_basin['iFlag_dam'] = 0
+        aConfig_basin['iFlag_disconnected'] = 0
+        aConfig_basin['lBasinID'] = i + 1
+        aConfig_basin['dLatitude_outlet_degree'] = -180
+        aConfig_basin['dLongitude_outlet_degree'] = 180
+        aConfig_basin['dAccumulation_threshold'] = -90
+        aConfig_basin['dThreshold_small_river'] = 90        
+        aConfig_basin['sFilename_dam'] = str(Path(sWorkspace_input_in)  /  'ICoM_dams.csv')        
+        aConfig_basin['sFilename_flowline_filter'] = str(Path(sWorkspace_input_in)  /  'streamord7above.shp')        
+        aConfig_basin['sFilename_flowline_raw'] = str(Path(sWorkspace_input_in)  /  'allflowline.shp')        
+        aConfig_basin['sFilename_flowline_topo'] = str(Path(sWorkspace_input_in)  /  'flowline.csv')        
+        aConfig_basin['sWorkspace_output_basin'] = str(Path(sWorkspace_output_in) / sBasin )
+        pBasin = pybasin(aConfig_basin)    
+        aBasin_out.append(pBasin)        
+        pass
+        
+    #export basin config to a file     
+    with open(sFilename_basins_json, 'w', encoding='utf-8') as f:
+        sJson = json.dumps([json.loads(ob.tojson()) for ob in aBasin_out], indent = 4)        
+        f.write(sJson)    
+        f.close()
+
+    return aBasin_out
+
+def pyflowline_create_template_configuration_file(
+        sFilename_json,
+        sWorkspace_input,
+        sWorkspace_output,
+        iFlag_standalone_in=None,
+        iFlag_use_mesh_dem_in=None,
+        iCase_index_in=None,
+        dResolution_degree_in=None,
+        dResolution_meter_in=None,
+        sDate_in=None,
+        sMesh_type_in=None,
+        sModel_in=None):
+    """generate pyflowline config template file
+
+    Args:
+        sFilename_json (str or Path): _description_.
+        sWorkspace_input (str or Path): _description_.
+        sWorkspace_output (str or Path): _description_.
+        iFlag_standalone_in (int, optional): _description_. Defaults to None.
+        iFlag_use_mesh_dem_in (int, optional): _description_. Defaults to None.        
+        iCase_index_in (int, optional): _description_. Defaults to None.
+        dResolution_degree_in (float, optional): _description_. Defaults to None.
+        dResolution_meter_in (float, optional): _description_. Defaults to None.
+        sDate_in (str, optional): _description_. Defaults to None.
+        sMesh_type_in (str, optional): _description_. Defaults to None.
+        sModel_in (str, optional): _description_. Defaults to None.
+
+    Returns:
+        _type_: _description_
+    """
+
+    # Ensure input pathnames are strings
+    sFilename_json = str(sFilename_json)
+    sWorkspace_input = str(sWorkspace_input)
+    sWorkspace_output = str(sWorkspace_output)
+
+    if os.path.exists(sFilename_json):         
+        os.remove(sFilename_json)
+
+    if iCase_index_in is not None:        
+        iCase_index = iCase_index_in
+    else:       
+        iCase_index = 1
+    
+    if iFlag_standalone_in is not None:        
+        iFlag_standalone = iFlag_standalone_in
+    else:       
+        iFlag_standalone = 1
+
+    if iFlag_use_mesh_dem_in is not None:        
+        iFlag_use_mesh_dem = iFlag_use_mesh_dem_in
+    else:       
+        iFlag_use_mesh_dem = 0
+
+    #if iFlag_use_shapefile_extent_in is not None:        
+    #    iFlag_use_shapefile_extent = iFlag_use_shapefile_extent_in
+    #else:       
+    #    iFlag_use_shapefile_extent = 0
+
+    if sMesh_type_in is not None:
+        sMesh_type = sMesh_type_in
+    else:
+        sMesh_type = 'hexagon'
+        pass
+    if sDate_in is not None:
+        sDate = sDate_in
+    else:
+        sDate = '20220202'
+        pass
+
+    nBasin = 1
+
+    # Use a dict to initialize the class
+    aConfig = {}
+
+    #aConfig['iFlag_use_shapefile_extent'] = iFlag_use_shapefile_extent 
+    aConfig['iFlag_use_mesh_dem'] = iFlag_use_mesh_dem
+    aConfig['iFlag_save_mesh'] = 1
+    aConfig['iFlag_simplification']=1
+    aConfig['iFlag_create_mesh']=1
+    aConfig['iFlag_intersect']=1
+    aConfig['iFlag_resample_method']=1
+    aConfig['iFlag_global']=0
+    aConfig['iFlag_multiple_outlet']=0
+    aConfig['iFlag_elevation_profile']=1
+    aConfig['iFlag_rotation']=0
+    aConfig['iFlag_stream_burning_topology']=1
+    aConfig['iFlag_save_elevation']=1     
+    aConfig['nOutlet'] = nBasin
+    aConfig['dResolution_degree'] = 0.5
+    aConfig['dResolution_meter'] = 50000
+    aConfig['dLongitude_left'] = -180
+    aConfig['dLongitude_right'] = 180
+    aConfig['dLatitude_bot'] = -90
+    aConfig['dLatitude_top'] = 90
+    aConfig['sFilename_model_configuration'] = sFilename_json 
+    aConfig['sWorkspace_input'] = sWorkspace_input    
+    aConfig['sWorkspace_output'] = sWorkspace_output        
+    aConfig['sRegion'] = 'susquehanna'
+    aConfig['sModel'] = 'pyflowline'
+    aConfig['iCase_index'] = iCase_index   
+    aConfig['sMesh_type'] = sMesh_type
+    aConfig['sJob'] = 'pyflowline'
+    aConfig['sDate']= sDate
+    #full path
+    aConfig['sFilename_mesh_netcdf'] = str(Path(sWorkspace_input)  /  'lnd_cull_mesh.nc')
+    aConfig['sFilename_spatial_reference'] =  str(Path(sWorkspace_input)  /  'boundary_proj.shp')   
+    #relative path
+    aConfig['flowline_info'] = 'flowline_info.json'
+    aConfig['sFilename_mesh_info'] = 'mesh_info.json'
+    aConfig['sFilename_elevation'] = 'elevation.json' 
+    oModel = flowlinecase(aConfig)
+
+    # Generate basin
+    sDirname = os.path.dirname(sFilename_json)
+    sFilename =  Path(sFilename_json).stem + '_basins.json'
+    sFilename_basins_json = os.path.join(sDirname, sFilename)   
+
+    aBasin = pyflowline_create_template_basin_configuration_file(
+        sFilename_basins_json,
+        nBasin,
+        sWorkspace_input,
+        oModel.sWorkspace_output)
+
+    oModel.aBasin = aBasin
+    oModel.sFilename_basins = sFilename_basins_json
+    oModel.export_config_to_json(sFilename_json)    
+
+    return oModel
diff --git a/pyflowline/configuration/path_manager.py b/pyflowline/configuration/path_manager.py
new file mode 100644
index 00000000..8318a6d7
--- /dev/null
+++ b/pyflowline/configuration/path_manager.py
@@ -0,0 +1,57 @@
+from pathlib import Path
+from pkg_resources import resource_filename
+"""
+path_manager
+All paths are created with Path and call resolve()
+"""
+
+def join_project_path(*args: str) -> Path:
+    """
+    Join components of a path relative to the project root.
+
+    Args:
+        *args: Variable number of string arguments representing path components.
+
+    Returns:
+        Path: Absolute path joined with the project root.
+    """
+    project_root = pyflowline_project_root()
+    return project_root.joinpath(*args)
+
+def pyflowline_project_root() -> Path:
+    """
+    Attempt to find the root path of the project, first by looking for a 'setup.py' file, and then using pkg_resources to find the installation path.
+    """
+    try:
+        return root_path_from_setup_file()
+    except FileNotFoundError:
+        return root_path_from_pyflowline_package_root()
+    except Exception as e:
+        raise RuntimeError(f"An unexpected error occurred while locating the project root: {str(e)}")
+
+def pyflowline_package_root(package_name='pyflowline') -> Path:
+    """Get the installation root directory of the package."""
+    return Path(resource_filename(package_name, '')).resolve()
+
+def root_path_from_setup_file() -> Path:
+    """Return the top-level (root) path of the pyflowline project.
+    This function navigates upwards from this file's path until it finds a directory
+    with a specific marker (e.g., a setup.py file) indicating the root of the project.
+    """
+    this_path = Path(__file__).resolve()
+    for parent in this_path.parents:
+        if (parent / 'setup.py').exists():
+            return parent
+    raise FileNotFoundError("setup.py not found. Is this the right project structure?")
+
+def root_path_from_pyflowline_package_root() -> Path:
+    """Get the installation root directory of the package."""
+    return pyflowline_package_root().parent
+
+# Example additional function for getting root from this module's location, adjusted for best practice
+def root_path_from_path_manager_location() -> Path:
+    """Return the top-level (root) path of the pyflowline project
+    This function assumes the module is located at a fixed level:
+        root_path/pyflowline/configuration/path_manager
+    """
+    return Path(__file__).resolve().parents[2]
diff --git a/pyflowline/configuration/pyflowline_create_template_configuration_file.py b/pyflowline/configuration/pyflowline_create_template_configuration_file.py
index 49d8baa1..81ad2bc7 100644
--- a/pyflowline/configuration/pyflowline_create_template_configuration_file.py
+++ b/pyflowline/configuration/pyflowline_create_template_configuration_file.py
@@ -6,19 +6,28 @@
 from pyflowline.classes.pycase import flowlinecase
 from pyflowline.classes.basin import pybasin
 
-def pyflowline_create_basin_template_configuration_file(sFilename_basins_json, nBasin, sWorkspace_input_in, sWorkspace_output_in):
+def pyflowline_create_template_basin_configuration_file(
+        sFilename_basins_json,
+        nBasin,
+        sWorkspace_input_in,
+        sWorkspace_output_in):
     """generate basin configuration
 
     Args:
-        sFilename_basins_json (str): the filename
+        sFilename_basins_json (str or Path): the filename
         nBasin (int): the total number of basin
-        sWorkspace_input_in (str): the input data path
-        sWorkspace_output (str): the output path
+        sWorkspace_input_in (str or Path): the input data path
+        sWorkspace_output (str or Path): the output path
 
     Returns:
         basin: a basin object
     """
 
+    # Ensure input pathnames are strings
+    sFilename_basins_json = str(sFilename_basins_json)
+    sWorkspace_input_in = str(sWorkspace_input_in)
+    sWorkspace_output_in = str(sWorkspace_output_in)
+
     aBasin_out = list()
     for i in range(nBasin):
         sBasin =  "{:03d}".format(i+1)   
@@ -47,14 +56,24 @@ def pyflowline_create_basin_template_configuration_file(sFilename_basins_json, n
 
     return aBasin_out
 
-def pyflowline_create_template_configuration_file(sFilename_json, \
-    sWorkspace_input, sWorkspace_output, iFlag_standalone_in=None, \
-        iFlag_use_mesh_dem_in=None,  iCase_index_in=None,  dResolution_degree_in = None, dResolution_meter_in = None,   sDate_in = None,  sMesh_type_in = None,    sModel_in = None):
+def pyflowline_create_template_configuration_file(
+        sFilename_json,
+        sWorkspace_input,
+        sWorkspace_output,
+        iFlag_standalone_in=None,
+        iFlag_use_mesh_dem_in=None,
+        iCase_index_in=None,
+        dResolution_degree_in=None,
+        dResolution_meter_in=None,
+        sDate_in=None,
+        sMesh_type_in=None,
+        sModel_in=None):
     """generate pyflowline config template file
 
     Args:
-        sFilename_json (str): _description_
-        sWorkspace_input (str): _description_
+        sFilename_json (str or Path): _description_.
+        sWorkspace_input (str or Path): _description_.
+        sWorkspace_output (str or Path): _description_.
         iFlag_standalone_in (int, optional): _description_. Defaults to None.
         iFlag_use_mesh_dem_in (int, optional): _description_. Defaults to None.        
         iCase_index_in (int, optional): _description_. Defaults to None.
@@ -63,11 +82,16 @@ def pyflowline_create_template_configuration_file(sFilename_json, \
         sDate_in (str, optional): _description_. Defaults to None.
         sMesh_type_in (str, optional): _description_. Defaults to None.
         sModel_in (str, optional): _description_. Defaults to None.
-        sWorkspace_output_in (str, optional): _description_. Defaults to None.
 
     Returns:
         _type_: _description_
     """
+
+    # Ensure input pathnames are strings
+    sFilename_json = str(sFilename_json)
+    sWorkspace_input = str(sWorkspace_input)
+    sWorkspace_output = str(sWorkspace_output)
+
     if os.path.exists(sFilename_json):         
         os.remove(sFilename_json)
 
@@ -101,15 +125,12 @@ def pyflowline_create_template_configuration_file(sFilename_json, \
     else:
         sDate = '20220202'
         pass
-    
-    
+
     nBasin = 1
 
- 
-    #use a dict to initialize the class
+    # Use a dict to initialize the class
     aConfig = {}
-    
-    
+
     #aConfig['iFlag_use_shapefile_extent'] = iFlag_use_shapefile_extent 
     aConfig['iFlag_use_mesh_dem'] = iFlag_use_mesh_dem
     aConfig['iFlag_save_mesh'] = 1
@@ -130,7 +151,7 @@ def pyflowline_create_template_configuration_file(sFilename_json, \
     aConfig['dLongitude_right'] = 180
     aConfig['dLatitude_bot'] = -90
     aConfig['dLatitude_top'] = 90
-    aConfig['sFilename_model_configuration']  = sFilename_json 
+    aConfig['sFilename_model_configuration'] = sFilename_json 
     aConfig['sWorkspace_input'] = sWorkspace_input    
     aConfig['sWorkspace_output'] = sWorkspace_output        
     aConfig['sRegion'] = 'susquehanna'
@@ -147,14 +168,20 @@ def pyflowline_create_template_configuration_file(sFilename_json, \
     aConfig['sFilename_mesh_info'] = 'mesh_info.json'
     aConfig['sFilename_elevation'] = 'elevation.json' 
     oModel = flowlinecase(aConfig)
-    #generate basin
+
+    # Generate basin
     sDirname = os.path.dirname(sFilename_json)
     sFilename =  Path(sFilename_json).stem + '_basins.json'
     sFilename_basins_json = os.path.join(sDirname, sFilename)   
-    aBasin = pyflowline_create_basin_template_configuration_file(sFilename_basins_json, nBasin, sWorkspace_input,  oModel.sWorkspace_output)
+
+    aBasin = pyflowline_create_template_basin_configuration_file(
+        sFilename_basins_json,
+        nBasin,
+        sWorkspace_input,
+        oModel.sWorkspace_output)
+
     oModel.aBasin = aBasin
     oModel.sFilename_basins = sFilename_basins_json
     oModel.export_config_to_json(sFilename_json)    
 
     return oModel
-
diff --git a/pyflowline/configuration/pyflowline_read_configuration_file.py b/pyflowline/configuration/pyflowline_read_configuration_file.py
index 1226714b..65b0ca7e 100644
--- a/pyflowline/configuration/pyflowline_read_configuration_file.py
+++ b/pyflowline/configuration/pyflowline_read_configuration_file.py
@@ -3,10 +3,26 @@
 import datetime
 import json
 from pyflowline.classes.pycase import flowlinecase
-pDate = datetime.datetime.today()
-sDate_default = "{:04d}".format(
-    pDate.year) + "{:02d}".format(pDate.month) + "{:02d}".format(pDate.day)
 
+pDate = datetime.datetime.today()
+sDate_default = "{:04d}".format(pDate.year) + \
+    "{:02d}".format(pDate.month) + "{:02d}".format(pDate.day)
+
+<<<<<<<< HEAD:pyflowline/configuration/read_configuration_file.py
+def pyflowline_read_configuration_file(
+        sFilename_configuration_in,
+        iFlag_standalone_in=None,
+        iFlag_use_mesh_dem_in=None,
+        iCase_index_in=None,
+        iResolution_index_in=None,
+        dResolution_degree_in=None,
+        dResolution_meter_in=None,
+        sMesh_type_in=None,
+        sModel_in=None,
+        sDate_in=None,
+        sDggrid_type_in=None,
+        sWorkspace_output_in=None):
+========
 
 def pyflowline_read_configuration_file(sFilename_configuration_in,
                                              iFlag_standalone_in=None,
@@ -20,10 +36,11 @@ def pyflowline_read_configuration_file(sFilename_configuration_in,
                                              sDate_in=None,
                                               sDggrid_type_in = None,
                                              sWorkspace_output_in=None):
+>>>>>>>> main:pyflowline/configuration/pyflowline_read_configuration_file.py
     """read a model configuration
 
     Args:
-        sFilename_configuration_in (str): _description_
+        sFilename_configuration_in (str or Path): _description_
         iFlag_standalone_in (int, optional): _description_. Defaults to None.
         iFlag_use_mesh_dem_in (int, optional): _description_. Defaults to None.
         iCase_index_in (int, optional): _description_. Defaults to None.
@@ -32,11 +49,16 @@ def pyflowline_read_configuration_file(sFilename_configuration_in,
         sMesh_type_in (str, optional): _description_. Defaults to None.
         sModel_in (str, optional): _description_. Defaults to None.
         sDate_in (str, optional): _description_. Defaults to None.
-        sWorkspace_output_in (str, optional): _description_. Defaults to None.
+        sWorkspace_output_in (str or Path, optional): _description_. Defaults to None.
 
     Returns:
         _type_: _description_
     """
+    # Ensure input filenames are strings
+    if isinstance(sFilename_configuration_in, Path):
+        sFilename_configuration_in = str(sFilename_configuration_in)
+    if isinstance(sWorkspace_output_in, Path):
+        sWorkspace_output_in = str(sWorkspace_output_in)
 
     if not os.path.isfile(sFilename_configuration_in):
         print(sFilename_configuration_in + ' does not exist')
@@ -116,8 +138,9 @@ def pyflowline_read_configuration_file(sFilename_configuration_in,
         # try to create this output folder first using
 
     try:
-        print(sWorkspace_output)
+        print("Creating the specified output workspace (if it does not exist): \n", sWorkspace_output)
         Path(sWorkspace_output).mkdir(parents=True, exist_ok=True)
+        print("The specified output workspace is: \n", sWorkspace_output)
     except ValueError:
         print("The specified output workspace cannot be created!")
         exit
diff --git a/pyflowline/configuration/read_configuration_file.py b/pyflowline/configuration/read_configuration_file.py
new file mode 100644
index 00000000..65b0ca7e
--- /dev/null
+++ b/pyflowline/configuration/read_configuration_file.py
@@ -0,0 +1,173 @@
+import os
+from pathlib import Path
+import datetime
+import json
+from pyflowline.classes.pycase import flowlinecase
+
+pDate = datetime.datetime.today()
+sDate_default = "{:04d}".format(pDate.year) + \
+    "{:02d}".format(pDate.month) + "{:02d}".format(pDate.day)
+
+<<<<<<<< HEAD:pyflowline/configuration/read_configuration_file.py
+def pyflowline_read_configuration_file(
+        sFilename_configuration_in,
+        iFlag_standalone_in=None,
+        iFlag_use_mesh_dem_in=None,
+        iCase_index_in=None,
+        iResolution_index_in=None,
+        dResolution_degree_in=None,
+        dResolution_meter_in=None,
+        sMesh_type_in=None,
+        sModel_in=None,
+        sDate_in=None,
+        sDggrid_type_in=None,
+        sWorkspace_output_in=None):
+========
+
+def pyflowline_read_configuration_file(sFilename_configuration_in,
+                                             iFlag_standalone_in=None,
+                                             iFlag_use_mesh_dem_in=None,
+                                             iCase_index_in=None,
+                                             iResolution_index_in = None,
+                                             dResolution_degree_in=None,
+                                             dResolution_meter_in=None,
+                                             sMesh_type_in=None,
+                                             sModel_in=None,
+                                             sDate_in=None,
+                                              sDggrid_type_in = None,
+                                             sWorkspace_output_in=None):
+>>>>>>>> main:pyflowline/configuration/pyflowline_read_configuration_file.py
+    """read a model configuration
+
+    Args:
+        sFilename_configuration_in (str or Path): _description_
+        iFlag_standalone_in (int, optional): _description_. Defaults to None.
+        iFlag_use_mesh_dem_in (int, optional): _description_. Defaults to None.
+        iCase_index_in (int, optional): _description_. Defaults to None.
+        dResolution_degree_in (float, optional): _description_. Defaults to None.
+        dResolution_meter_in (float, optional): _description_. Defaults to None.
+        sMesh_type_in (str, optional): _description_. Defaults to None.
+        sModel_in (str, optional): _description_. Defaults to None.
+        sDate_in (str, optional): _description_. Defaults to None.
+        sWorkspace_output_in (str or Path, optional): _description_. Defaults to None.
+
+    Returns:
+        _type_: _description_
+    """
+    # Ensure input filenames are strings
+    if isinstance(sFilename_configuration_in, Path):
+        sFilename_configuration_in = str(sFilename_configuration_in)
+    if isinstance(sWorkspace_output_in, Path):
+        sWorkspace_output_in = str(sWorkspace_output_in)
+
+    if not os.path.isfile(sFilename_configuration_in):
+        print(sFilename_configuration_in + ' does not exist')
+        return
+
+    # Opening JSON file
+    with open(sFilename_configuration_in) as json_file:
+        aConfig = json.load(json_file)
+
+    if iCase_index_in is not None:
+        iCase_index = iCase_index_in
+    else:
+        iCase_index = int(aConfig['iCase_index'])
+
+    if iResolution_index_in is not None:
+        iResolution_index = iResolution_index_in
+    else:
+        if "iResolution_index" in aConfig:
+            iResolution_index =  int( aConfig['iResolution_index'])
+        else:
+            iResolution_index = 10
+
+        pass
+
+    if iFlag_standalone_in is not None:
+        iFlag_standalone = iFlag_standalone_in
+    else:
+        iFlag_standalone = int(aConfig['iFlag_standalone'])
+
+    if iFlag_use_mesh_dem_in is not None:
+        iFlag_use_mesh_dem = iFlag_use_mesh_dem_in
+    else:
+        iFlag_use_mesh_dem = int(aConfig['iFlag_use_mesh_dem'])
+
+    if sMesh_type_in is not None:
+        sMesh_type = sMesh_type_in
+    else:
+        sMesh_type = aConfig['sMesh_type']
+        pass
+
+    if sDggrid_type_in is not None:
+        sDggrid_type = sDggrid_type_in
+    else:
+        if "sDggrid_type" in aConfig:
+            sDggrid_type = aConfig["sDggrid_type"]
+        else:
+            sDggrid_type = 'ISEA3H'
+
+    if sModel_in is not None:
+        sModel = sModel_in
+    else:
+        sModel = aConfig['sModel']
+        pass
+
+    if sDate_in is not None:
+        sDate = sDate_in
+    else:
+        sDate = aConfig['sDate']
+        pass
+
+    if dResolution_degree_in is not None:
+        dResolution_degree = dResolution_degree_in
+    else:
+        dResolution_degree = aConfig['dResolution_degree']
+        pass
+
+    if dResolution_meter_in is not None:
+        dResolution_meter = dResolution_meter_in
+    else:
+        dResolution_meter = aConfig['dResolution_meter']
+        pass
+
+    if sWorkspace_output_in is not None:
+        sWorkspace_output = sWorkspace_output_in
+    else:
+        sWorkspace_output = aConfig['sWorkspace_output']
+        # try to create this output folder first using
+
+    try:
+        print("Creating the specified output workspace (if it does not exist): \n", sWorkspace_output)
+        Path(sWorkspace_output).mkdir(parents=True, exist_ok=True)
+        print("The specified output workspace is: \n", sWorkspace_output)
+    except ValueError:
+        print("The specified output workspace cannot be created!")
+        exit
+
+    aConfig['iCase_index'] = iCase_index
+    aConfig['iFlag_standalone'] = iFlag_standalone
+    aConfig['iFlag_use_mesh_dem'] = iFlag_use_mesh_dem
+    aConfig["iResolution_index"] = iResolution_index
+    aConfig['dResolution_degree'] = dResolution_degree
+    aConfig['dResolution_meter'] = dResolution_meter
+
+    aConfig['sDate'] = sDate
+    aConfig['sModel'] = sModel
+    aConfig['sMesh_type'] = sMesh_type
+    aConfig["sDggrid_type"] = sDggrid_type
+    aConfig['sWorkspace_output'] = sWorkspace_output
+
+    aConfig["sFilename_model_configuration"] = sFilename_configuration_in
+
+    # based on global variable, a few variables are calculate once
+    # calculate the modflow simulation period
+    # https://docs.python.org/3/library/datetime.html#datetime-objects
+
+    # aConfig
+
+    # simulation
+
+    oPyflowline = flowlinecase(aConfig)
+
+    return oPyflowline
diff --git a/tests/TestPathManager.py b/tests/TestPathManager.py
new file mode 100644
index 00000000..69757cd1
--- /dev/null
+++ b/tests/TestPathManager.py
@@ -0,0 +1,32 @@
+import unittest
+from pathlib import Path
+from pyflowline.configuration import path_manager
+
+class TestPathManager(unittest.TestCase):
+
+    def test_root_path_from_pyflowline_project_root(self):
+        """Test the dynamic root path identification using both setup.py and pkg_resources."""
+        expected_root = Path(__file__).resolve().parents[1]
+        returned_root = path_manager.pyflowline_project_root()
+        self.assertEqual(returned_root, expected_root, "Root path is incorrectly identified from pyflowline project root.")
+
+    def test_root_path_from_pyflowline_package_root(self):
+        """Test the package root path identification using pkg_resources."""
+        expected_root = Path(__file__).resolve().parents[1]
+        returned_root = path_manager.root_path_from_pyflowline_package_root()
+        self.assertEqual(returned_root, expected_root, "Root path is incorrectly identified from pyflowline package root.")
+
+    def test_root_path_from_setup_file(self):
+        """Test the root path identification specifically by locating setup.py."""
+        expected_root = Path(__file__).resolve().parents[1]
+        returned_root = path_manager.root_path_from_setup_file()
+        self.assertEqual(returned_root, expected_root, "Root path is incorrectly identified from setup file location.")
+
+    def test_root_path_from_path_manager_location(self):
+        """Test the root path based on the relative location of the path_manager module."""
+        expected_root = Path(__file__).resolve().parents[1]
+        returned_root = path_manager.root_path_from_path_manager_location()
+        self.assertEqual(returned_root, expected_root, "Root path is incorrectly identified from path_manager location.")
+
+if __name__ == '__main__':
+    unittest.main()
diff --git a/tests/test_installation.py b/tests/test_installation.py
index efc15401..e982f459 100644
--- a/tests/test_installation.py
+++ b/tests/test_installation.py
@@ -9,4 +9,11 @@ def test_installation():
 
 test_installation()
 
-
+# This might be useful:
+# import pkg_resources
+# def test_pkg_installation():
+#     try:
+#         distribution = pkg_resources.get_distribution('pyflowline')
+#         print("Pyflowline is installed at:", distribution.location)
+#     except pkg_resources.DistributionNotFound:
+#         print("Pyflowline is not installed.")