diff --git a/common/Rias_UWWTP.xlsx b/common/Rias_UWWTP.xlsx
new file mode 100644
index 0000000..d820225
Binary files /dev/null and b/common/Rias_UWWTP.xlsx differ
diff --git a/common/__init__.py b/common/__init__.py
new file mode 100644
index 0000000..3f5aad4
--- /dev/null
+++ b/common/__init__.py
@@ -0,0 +1,3 @@
+from .inout import read_input
+from .boundarybox import BoundaryBox
+from .vector import Vector, NumpyVector
diff --git a/common/__pycache__/__init__.cpython-39.pyc b/common/__pycache__/__init__.cpython-39.pyc
new file mode 100644
index 0000000..49d9a5c
Binary files /dev/null and b/common/__pycache__/__init__.cpython-39.pyc differ
diff --git a/common/__pycache__/boundarybox.cpython-39.pyc b/common/__pycache__/boundarybox.cpython-39.pyc
new file mode 100644
index 0000000..fec69f3
Binary files /dev/null and b/common/__pycache__/boundarybox.cpython-39.pyc differ
diff --git a/common/__pycache__/drawcurrents.cpython-39.pyc b/common/__pycache__/drawcurrents.cpython-39.pyc
new file mode 100644
index 0000000..498ba01
Binary files /dev/null and b/common/__pycache__/drawcurrents.cpython-39.pyc differ
diff --git a/common/__pycache__/drawmap.cpython-39.pyc b/common/__pycache__/drawmap.cpython-39.pyc
new file mode 100644
index 0000000..7cf79b1
Binary files /dev/null and b/common/__pycache__/drawmap.cpython-39.pyc differ
diff --git a/common/__pycache__/inout.cpython-39.pyc b/common/__pycache__/inout.cpython-39.pyc
new file mode 100644
index 0000000..3f2b838
Binary files /dev/null and b/common/__pycache__/inout.cpython-39.pyc differ
diff --git a/common/__pycache__/map.cpython-39.pyc b/common/__pycache__/map.cpython-39.pyc
new file mode 100644
index 0000000..494d34f
Binary files /dev/null and b/common/__pycache__/map.cpython-39.pyc differ
diff --git a/common/__pycache__/vector.cpython-39.pyc b/common/__pycache__/vector.cpython-39.pyc
new file mode 100644
index 0000000..3b08a8a
Binary files /dev/null and b/common/__pycache__/vector.cpython-39.pyc differ
diff --git a/common/boundarybox.py b/common/boundarybox.py
new file mode 100644
index 0000000..e8c2824
--- /dev/null
+++ b/common/boundarybox.py
@@ -0,0 +1,20 @@
+class BoundaryBox:
+    lat_min: float
+    lat_max: float
+    lon_min: float
+    lon_max: float
+
+    def __init__(self, lat_min=-90., lat_max=90., lon_min=-180., lon_max=180.):
+        self.lat_min = lat_min
+        self.lat_max = lat_max
+        self.lon_min = lon_min
+        self.lon_max = lon_max
+
+    def middle_lat(self):
+        x = self.lat_min + 0.5 * (self.lat_max - self.lat_min)
+        return x
+
+    def middle_lon(self):
+        return self.lon_min + 0.5 * (self.lon_max - self.lon_min)
+
+
diff --git a/common/drawcurrents.py b/common/drawcurrents.py
new file mode 100644
index 0000000..04ab4ae
--- /dev/null
+++ b/common/drawcurrents.py
@@ -0,0 +1,81 @@
+
+
+def drawcurrents(coordinates_rank, nx, ny, scale, resolution,
+                 level, time, lat, lon, ust, vst, mod,
+                 file_out, title, style, boundary_box):
+    """
+
+    :return:
+    """
+    import numpy as np
+    import matplotlib.pyplot as plt
+    from mpl_toolkits.basemap import Basemap
+
+    fig = plt.figure(figsize=(10, 10))
+    ax = fig.add_subplot(111)
+
+    dx = 0.01
+    middle_lon = boundary_box.middle_lon()
+    middle_lat = boundary_box.middle_lat()
+    m = Basemap(llcrnrlon=boundary_box.lon_min - dx,
+                llcrnrlat=boundary_box.lat_min - dx,
+                urcrnrlon=boundary_box.lon_max + dx,
+                urcrnrlat=boundary_box.lat_max + dx,
+                resolution=resolution, projection='tmerc', lon_0=middle_lon, lat_0=middle_lat)
+
+    m.drawcoastlines()
+    m.fillcontinents(color='grey', lake_color='aqua')
+
+
+    #m.drawmapboundary(fill_color='aqua')
+
+    if coordinates_rank == 1:
+        lon, lat = np.meshgrid(lon, lat)
+
+    x, y = m(lon, lat)
+
+    # draw filled contours.
+    clevs = [0, 0.1, 0.2, 0.3, 0.4, 0.5]
+    # clevs = [0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1]
+    if style == 'contour':
+        cs = m.contourf(x, y, mod, clevs, cmap=plt.cm.jet)
+        plt.quiver(x[::nx, ::ny], y[::nx, ::ny], ust[::nx, ::ny], vst[::nx, ::ny], scale=scale)
+
+    if style == 'cvector':
+        clim = [0., 0.5]
+        cs = plt.quiver(x[::nx, ::ny], y[::nx, ::ny], ust[::nx, ::ny], vst[::nx, ::ny], mod[::nx, ::ny], clim=clim, scale=scale)
+
+    if style == 'windbarbs':
+        clim = [0., 10]
+        cs = plt.barbs(x[::nx, ::ny], y[::nx, ::ny], ust[::nx, ::ny], vst[::nx, ::ny], mod[::nx, ::ny], clim=clim,cmap=plt.cm.jet)
+
+    if style == 'stream':
+        print(lons.shape)
+        print(x[::, 1])
+        print(ust.shape)
+        xx = x[1::, ::]
+        xxx = xx[0]
+        yyy = y[::, 1::][1]
+
+        u = ust[::, ::]
+        v = vst[::, ::]
+        print(len(xxx))
+        print(len(yyy))
+        print(len(u))
+        print(len(v))
+
+        cs = plt.streamplot(lon, lat, u, v)
+
+    cbar = m.colorbar(cs, location='bottom', pad="5%")
+    cbar.ax.tick_params(labelsize=9)
+    cbar.set_label('m/s')
+
+    # add title
+    plt.title(title)
+
+    fig.savefig(file_out, dpi=100, facecolor='w', edgecolor='w', format='png',
+                transparent=False,  bbox_inches='tight', pad_inches=0.1)
+    plt.clf()
+    plt.close(fig)
+
+    return
\ No newline at end of file
diff --git a/common/drawmap.json b/common/drawmap.json
new file mode 100644
index 0000000..58559c1
--- /dev/null
+++ b/common/drawmap.json
@@ -0,0 +1,28 @@
+{
+    "path_in": "",
+    "path_out": "",
+    "file_in": "",
+    "file_out": "",
+    "nx": 1,
+    "ny": 1,
+    "vector": 0,
+    "scalar": 1,
+    "u": "eastward_velocity",
+    "v": "northward_velocity",
+    "scalar_magnitude": "",
+    "resolution": "i",
+    "scale": 5,
+    "n_time": "all",
+    "n_level": 9,
+    "title": "",
+    "style": "cvector",
+    "limits": [
+        42.08,
+        42.45,
+        -9.004,
+        -8.601
+    ],
+    "wms_url": "https://ows.emodnet-bathymetry.eu/wms?",
+    "wms_layers": ["emodnet:mean_atlas_land", "coastlines"]
+
+}
\ No newline at end of file
diff --git a/common/drawmap.py b/common/drawmap.py
new file mode 100644
index 0000000..36b1cd7
--- /dev/null
+++ b/common/drawmap.py
@@ -0,0 +1,252 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+
+# Para el tratamiento de datos
+import numpy as np
+import pandas as pd
+# Para las representaciones
+import matplotlib.pyplot as plt
+from matplotlib.cm import jet
+# Para los directorios
+import os
+# Importamos funciones de otros programas
+from common.readers.reader_factory import read_factory
+from common import read_input
+from common.boundarybox import BoundaryBox
+from common.drawcurrents import drawcurrents
+from common.map import Map
+
+def read_inputs(input_file):
+    """Read keywords for options"""
+    input_keys = ['path_in',
+                  'file_in',
+                  'path_out'
+                  'file_out',
+                  'nx',
+                  'ny',
+                  'resolution',
+                  'scale',
+                  'n_time',
+                  'n_level',
+                  'title',
+                  'style',
+                  'limits',
+                  'vector',
+                  'scalar',
+                  'wms_url',
+                  'wms_layers']
+    return read_input(input_file, input_keys)
+
+
+def draw_map_vector(inputs, n):
+    """draw 1 maps of a day"""
+    draw_map = Vector(inputs)
+    draw_map.read_head()
+    draw_map.create_title(n)
+    draw_map.reader_by_time()
+    draw_map.draw()
+
+
+def draw_map_scalar(inputs, n):
+    draw_map = Scalar(inputs)
+    draw_map.read_head()
+    draw_map.create_title(n)
+
+    if n is None:
+        draw_map.reader_no_time()
+    else:
+        draw_map.reader_by_time()
+    draw_map.draw()
+
+
+def draw_map_24(inputs,start,end,ref_lon, ref_lat, key, draw_scale):
+    """draw 24+1 maps of a day"""
+
+    if inputs['scalar']:
+        draw_map = Scalar(inputs)
+    elif inputs['vector']:
+        draw_map = Vector(inputs)
+    
+    draw_map.read_head()
+    draw_map.options.time = 0 # creo que esto sobra
+
+    for n in range(draw_map.reader.ini_ntime+start, draw_map.reader.ini_ntime+end):
+        draw_map.create_title(n)
+        draw_map.options.time = n
+        draw_map.reader_by_time()
+        draw_map.draw(ref_lon, ref_lat, key, draw_scale)
+
+
+class OptionsMap:
+
+    def __init__(self, inputs):
+
+        self.file_path_in = inputs['path_in']
+        self.file_path_out = inputs['path_out']
+        self.file_in = inputs['file_in']
+        self.file_name = os.path.join(self.file_path_in, self.file_in)
+        self.file_hdf_out = inputs['file_out']
+        self.file_out = os.path.join(self.file_path_out, self.file_hdf_out)
+
+        self.nx = inputs['nx']
+        self.ny = inputs['ny']
+        self.scale = inputs['scale']
+        self.resolution = inputs['resolution']
+        self.style = inputs['style']
+        self.title = inputs['title']
+
+        self.level = inputs['n_level']
+        self.time = inputs['n_time']
+        limits = inputs['limits']
+        self.boundary_box = BoundaryBox(limits[0], limits[1], limits[2], limits[3])
+
+        self.vector_bool = inputs['vector']
+        self.scalar_bool = inputs['scalar']
+
+        self.wms_url = inputs['wms_url']
+        self.wms_layers = inputs['wms_layers']
+
+
+class DrawMap:
+
+    def __init__(self, inputs):
+        self.options = OptionsMap(inputs)
+        self.reader = self.get_reader(self.options.file_name)
+
+    def read_head(self):
+
+        with self.reader.open():
+
+            lat = self.reader.latitudes
+            lon = self.reader.longitudes
+
+            if self.reader.coordinates_rank == 1:
+                self.lats = lat[0:self.reader.n_latitudes - 2]  # ATENCIóN:Esto y lo de abajo era -1, revisar
+                self.lons = lon[0:self.reader.n_longitudes - 2]
+            elif self.reader.coordinates_rank == 2:
+                self.lats = lat[0:self.reader.n_longitudes - 2, 0:self.reader.n_latitudes - 2]
+                self.lons = lon[0:self.reader.n_longitudes - 2, 0:self.reader.n_latitudes - 2]
+
+    def create_title(self, n_time):
+        if n_time is None:
+            n_time = 0
+        with self.reader.open():
+            data = self.reader.get_date(n_time)
+            data_str = data.strftime("%Y-%m-%d %H:%M UTC")
+            data_comp = data.strftime("%Y%m%d%H%M")
+            self.title_full = self.options.title + " " + data_str
+            self.file_out_full = self.options.file_out + '_' + data_comp + '.png'
+
+    def get_reader(self, file_in):
+        print('Opening: {0}'.format(file_in))
+        factory = read_factory(file_in)
+        return factory.get_reader()
+
+
+class Vector(DrawMap):
+
+    def __init__(self, inputs):
+
+        super().__init__(inputs)
+
+        self.u_name = inputs['u']
+        self.v_name = inputs['v']
+        self.us = None
+        self.vs = None
+        self.modules = None
+
+    def reader_by_time(self):
+
+        with self.reader.open():
+            u = self.reader.get_variable(self.u_name, self.options.time)
+            v = self.reader.get_variable(self.v_name, self.options.time)
+
+            if len(u.shape) == 3:
+                self.us = u[self.options.level, :-1, :- 1]
+                self.vs = v[self.options.level, :-1, :-1]
+
+            elif len(u.shape) == 2:
+                self.us = u[:-1, :-1]
+                self.vs = v[:-1, :-1]
+
+            self.modules = pow((pow(self.us, 2) + pow(self.vs, 2)), .5)
+
+    def draw(self):
+        drawcurrents(self.reader.coordinates_rank, self.options.nx, self.options.ny,
+                     self.options.scale, self.options.resolution,
+                     self.options.level, self.options.time,
+                     self.lats, self.lons, self.us, self.vs, self.modules,
+                     self.file_out_full, self.title_full, self.options.style,
+                     self.options.boundary_box)
+
+
+class Scalar(DrawMap):
+
+    def __init__(self, inputs):
+
+        super().__init__(inputs)
+
+        self.scalar_name = inputs['scalar_magnitude']
+        self.scalars = None
+
+    def reader_by_time(self):
+
+        with self.reader.open():
+            scalar = self.reader.get_variable(self.scalar_name, self.options.time)
+
+            if len(scalar.shape) == 3:
+                self.scalars = scalar[self.options.level, :-1, :- 1]
+
+            elif len(scalar.shape) == 2:
+                self.scalars = scalar[:-1, :-1]
+
+    def reader_no_time(self):
+        with self.reader.open():
+            scalar = self.reader.get_var(self.scalar_name)
+
+            if len(scalar.shape) == 3:
+                self.scalars = scalar[self.options.level, :-1, :- 1]
+
+            elif len(scalar.shape) == 2:
+                self.scalars = scalar[:-1, :-1]
+
+    def draw(self,ref_lon, ref_lat, key, draw_scale):
+        """
+        :return:
+        """
+        # a Map object is created to deal with
+        mapa = Map()
+        mapa.set_extension(self.options.boundary_box.lat_min, self.options.boundary_box.lat_max,
+                           self.options.boundary_box.lon_min, self.options.boundary_box.lon_max)
+        fig = mapa.add_figure()
+        mapa.add_wms(self.options.wms_url, self.options.wms_layers)
+        #mapa.add_wms("https://wms.mapama.gob.es/sig/Agua/PVERT/2021/wms.aspx?", 
+                       #["Puntos de vertido de depuradoras urbanas (Q2021. Dir 91/271/CEE)"])
+        
+
+        #  Draw the scalars
+        if self.reader.coordinates_rank == 1:
+            longitudes, latitudes = np.meshgrid(self.lons, self.lats)
+        else:
+            longitudes, latitudes = self.lons, self.lats
+
+        im = mapa.ax.pcolormesh(longitudes, latitudes, self.scalars[:-1, :-1], shading='auto')  # since lats and lons are n-1xn-1 matrix
+        if draw_scale=='jet':
+            im.set_cmap(jet)
+        
+        # add title
+        mapa.ax.set_title(self.title_full)
+        # add point
+        path_discharge = os.getcwd() +'\\common\\Rias_UWWTP.xlsx'
+        discharge_lats = (pd.read_excel(path_discharge, sheet_name = key[1])['Lat'].values)
+        discharge_lons = (pd.read_excel(path_discharge, sheet_name = key[1])['Lon'].values)
+        mapa.ax.scatter(discharge_lons, discharge_lats, marker='.', s=150, color='darkorange',
+                        edgecolors = "black", linewidths = 1.)
+        mapa.ax.scatter(ref_lon, ref_lat, marker='*', s=100, color='red')
+        fig.colorbar(im, ax=mapa.ax)
+        #
+        fig.savefig(self.file_out_full, dpi=100, facecolor='w', edgecolor='w', format='png',
+                    transparent=False, bbox_inches='tight', pad_inches=0.1)
+        plt.clf()
+        plt.close(fig)
+        return
diff --git a/common/drawmap_copia.py b/common/drawmap_copia.py
new file mode 100644
index 0000000..05eacbc
--- /dev/null
+++ b/common/drawmap_copia.py
@@ -0,0 +1,252 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+
+# Para el tratamiento de datos
+import numpy as np
+import pandas as pd
+# Para las representaciones
+import matplotlib.pyplot as plt
+from matplotlib.cm import jet
+# Para los directorios
+import os
+# Importamos funciones de otros programas
+from common.readers.reader_factory import read_factory
+from common import read_input
+from common.boundarybox import BoundaryBox
+from common.drawcurrents import drawcurrents
+from common.map import Map
+
+def read_inputs(input_file):
+    """Read keywords for options"""
+    input_keys = ['path_in',
+                  'file_in',
+                  'path_out'
+                  'file_out',
+                  'nx',
+                  'ny',
+                  'resolution',
+                  'scale',
+                  'n_time',
+                  'n_level',
+                  'title',
+                  'style',
+                  'limits',
+                  'vector',
+                  'scalar',
+                  'wms_url',
+                  'wms_layers']
+    return read_input(input_file, input_keys)
+
+
+def draw_map_vector(inputs, n):
+    """draw 1 maps of a day"""
+    draw_map = Vector(inputs)
+    draw_map.read_head()
+    draw_map.create_title(n)
+    draw_map.reader_by_time()
+    draw_map.draw()
+
+
+def draw_map_scalar(inputs, n):
+    draw_map = Scalar(inputs)
+    draw_map.read_head()
+    draw_map.create_title(n)
+
+    if n is None:
+        draw_map.reader_no_time()
+    else:
+        draw_map.reader_by_time()
+    draw_map.draw()
+
+
+def draw_map_24(inputs,start,end,ref_lon, ref_lat):
+    """draw 24+1 maps of a day"""
+
+    if inputs['scalar']:
+        draw_map = Scalar(inputs)
+    elif inputs['vector']:
+        draw_map = Vector(inputs)
+    
+    draw_map.read_head()
+    draw_map.options.time = 0 # creo que esto sobra
+
+    for n in range(draw_map.reader.ini_ntime+start, draw_map.reader.ini_ntime+end):
+        draw_map.create_title(n)
+        draw_map.options.time = n
+        draw_map.reader_by_time()
+        draw_map.draw(ref_lon, ref_lat)
+
+
+class OptionsMap:
+
+    def __init__(self, inputs):
+
+        self.file_path_in = inputs['path_in']
+        self.file_path_out = inputs['path_out']
+        self.file_in = inputs['file_in']
+        self.file_name = os.path.join(self.file_path_in, self.file_in)
+        self.file_hdf_out = inputs['file_out']
+        self.file_out = os.path.join(self.file_path_out, self.file_hdf_out)
+
+        self.nx = inputs['nx']
+        self.ny = inputs['ny']
+        self.scale = inputs['scale']
+        self.resolution = inputs['resolution']
+        self.style = inputs['style']
+        self.title = inputs['title']
+
+        self.level = inputs['n_level']
+        self.time = inputs['n_time']
+        limits = inputs['limits']
+        self.boundary_box = BoundaryBox(limits[0], limits[1], limits[2], limits[3])
+
+        self.vector_bool = inputs['vector']
+        self.scalar_bool = inputs['scalar']
+
+        self.wms_url = inputs['wms_url']
+        self.wms_layers = inputs['wms_layers']
+
+
+class DrawMap:
+
+    def __init__(self, inputs):
+        self.options = OptionsMap(inputs)
+        self.reader = self.get_reader(self.options.file_name)
+
+    def read_head(self):
+
+        with self.reader.open():
+
+            lat = self.reader.latitudes
+            lon = self.reader.longitudes
+
+            if self.reader.coordinates_rank == 1:
+                self.lats = lat[0:self.reader.n_latitudes - 2]  # ATENCIóN:Esto y lo de abajo era -1, revisar
+                self.lons = lon[0:self.reader.n_longitudes - 2]
+            elif self.reader.coordinates_rank == 2:
+                self.lats = lat[0:self.reader.n_longitudes - 2, 0:self.reader.n_latitudes - 2]
+                self.lons = lon[0:self.reader.n_longitudes - 2, 0:self.reader.n_latitudes - 2]
+
+    def create_title(self, n_time):
+        if n_time is None:
+            n_time = 0
+        with self.reader.open():
+            data = self.reader.get_date(n_time)
+            data_str = data.strftime("%Y-%m-%d %H:%M UTC")
+            data_comp = data.strftime("%Y%m%d%H%M")
+            self.title_full = self.options.title + " " + data_str
+            self.file_out_full = self.options.file_out + '_' + data_comp + '.png'
+
+    def get_reader(self, file_in):
+        print('Opening: {0}'.format(file_in))
+        factory = read_factory(file_in)
+        return factory.get_reader()
+
+
+class Vector(DrawMap):
+
+    def __init__(self, inputs):
+
+        super().__init__(inputs)
+
+        self.u_name = inputs['u']
+        self.v_name = inputs['v']
+        self.us = None
+        self.vs = None
+        self.modules = None
+
+    def reader_by_time(self):
+
+        with self.reader.open():
+            u = self.reader.get_variable(self.u_name, self.options.time)
+            v = self.reader.get_variable(self.v_name, self.options.time)
+
+            if len(u.shape) == 3:
+                self.us = u[self.options.level, :-1, :- 1]
+                self.vs = v[self.options.level, :-1, :-1]
+
+            elif len(u.shape) == 2:
+                self.us = u[:-1, :-1]
+                self.vs = v[:-1, :-1]
+
+            self.modules = pow((pow(self.us, 2) + pow(self.vs, 2)), .5)
+
+    def draw(self):
+        drawcurrents(self.reader.coordinates_rank, self.options.nx, self.options.ny,
+                     self.options.scale, self.options.resolution,
+                     self.options.level, self.options.time,
+                     self.lats, self.lons, self.us, self.vs, self.modules,
+                     self.file_out_full, self.title_full, self.options.style,
+                     self.options.boundary_box)
+
+
+class Scalar(DrawMap):
+
+    def __init__(self, inputs):
+
+        super().__init__(inputs)
+
+        self.scalar_name = inputs['scalar_magnitude']
+        self.scalars = None
+
+    def reader_by_time(self):
+
+        with self.reader.open():
+            scalar = self.reader.get_variable(self.scalar_name, self.options.time)
+
+            if len(scalar.shape) == 3:
+                self.scalars = scalar[self.options.level, :-1, :- 1]
+
+            elif len(scalar.shape) == 2:
+                self.scalars = scalar[:-1, :-1]
+
+    def reader_no_time(self):
+        with self.reader.open():
+            scalar = self.reader.get_var(self.scalar_name)
+
+            if len(scalar.shape) == 3:
+                self.scalars = scalar[self.options.level, :-1, :- 1]
+
+            elif len(scalar.shape) == 2:
+                self.scalars = scalar[:-1, :-1]
+
+    def draw(self,ref_lon, ref_lat):
+        """
+        :return:
+        """
+        # a Map object is created to deal with
+        mapa = Map()
+        mapa.set_extension(self.options.boundary_box.lat_min, self.options.boundary_box.lat_max,
+                           self.options.boundary_box.lon_min, self.options.boundary_box.lon_max)
+        fig = mapa.add_figure()
+        mapa.add_wms(self.options.wms_url, self.options.wms_layers)
+        #mapa.add_wms("https://wms.mapama.gob.es/sig/Agua/PVERT/2021/wms.aspx?", 
+                       #["Puntos de vertido de depuradoras urbanas (Q2021. Dir 91/271/CEE)"])
+        
+
+        #  Draw the scalars
+        if self.reader.coordinates_rank == 1:
+            longitudes, latitudes = np.meshgrid(self.lons, self.lats)
+        else:
+            longitudes, latitudes = self.lons, self.lats
+
+        im = mapa.ax.pcolormesh(longitudes, latitudes, self.scalars[:-1, :-1], shading='auto')  # since lats and lons are n-1xn-1 matrix
+        #im.set_cmap(jet)
+        # add title
+        mapa.ax.set_title(self.title_full)
+        # add point
+        path_discharge = os.getcwd() +'\\common\\Rias_UWWTP.xlsx'
+
+        key = ['arousa', 'Arousa']
+        discharge_lats = (pd.read_excel(path_discharge, sheet_name = key[1])['Lat'].values)
+        discharge_lons = (pd.read_excel(path_discharge, sheet_name = key[1])['Lon'].values)
+        mapa.ax.scatter(discharge_lons, discharge_lats, marker='.', s=180, color='darkorange',
+                        edgecolors = "black", linewidths = 1.)
+        mapa.ax.scatter(ref_lon, ref_lat, marker='*', s=100, color='red')
+        fig.colorbar(im, ax=mapa.ax)
+        #
+        fig.savefig(self.file_out_full, dpi=100, facecolor='w', edgecolor='w', format='png',
+                    transparent=False, bbox_inches='tight', pad_inches=0.1)
+        plt.clf()
+        plt.close(fig)
+        return
diff --git a/common/hdf5_to_nc4.py b/common/hdf5_to_nc4.py
new file mode 100644
index 0000000..ac7ea3e
--- /dev/null
+++ b/common/hdf5_to_nc4.py
@@ -0,0 +1,164 @@
+# -*- coding: utf-8 -*-
+"""
+Created on Tue Mar 28 12:31:18 2023
+
+@author: josmo
+"""
+
+
+import h5py
+import xarray as xr
+import numpy as np
+import urllib
+import pandas as pd
+# Para las fechas
+from datetime import datetime,timedelta
+import os
+from tqdm import tqdm
+import glob
+
+# Rutas no hdf5
+
+# Abre o ficheiro lagranxiano e le
+
+
+def MOHID_hdf5_to_nc(file_in):
+   
+    rootGroup = 'Results/'
+    nameTimesGroup = '/Time/'
+    nameTimesNameGroup = '/Time/Time_'
+   
+    f = h5py.File(file_in, 'r')
+   
+    # Le as latitudes e lonxitudes e busca os nodos do cadrado
+   
+    latIn = f['/Grid/Latitude']
+    lonIn = f['/Grid/Longitude']
+   
+    rank = len(latIn.shape)
+   
+    if rank == 2:
+        lat = latIn[0,]
+        lon = lonIn[:, 0]
+   
+    if rank == 1:
+        lat = latIn
+        lon = lonIn
+   
+    # le os tempos pois o bucle vai ser por aqui
+    timesGroup = f[nameTimesGroup]
+    tempo = []
+    u_t = []
+    v_t = []
+   
+    for nameTime in timesGroup:
+       
+        numName = nameTime.split('_')[1]
+        valNumName = int(numName)
+        rootNameTime = nameTimesGroup + nameTime
+        time = f[rootNameTime]
+       
+        dataIn = datetime(int(time[0]), int(time[1]), int(time[2]), int(time[3]), int(time[4]), int(time[5]))
+        # nameData = '%04d/%02d/%02d %02d:%02d:%02d' % (int(time[0]),int(time[1]),int(time[2]),int(time[3]),int(time[4]),int(time[5]))
+       
+        tempo.append(dataIn)
+        nameU = rootGroup + "velocity U/velocity U_" + numName
+        nameV = rootGroup + "velocity V/velocity V_" + numName
+       
+        u = f[nameU]
+        v = f[nameV]
+        u_t.append(u[-1,:,:])
+        v_t.append(v[-1,:,:])
+   
+   
+    lon_c = (lon[1:]+lon[0:-1])/2.
+    lat_c = (lat[1:]+lat[0:-1])/2.
+    U = np.stack(u_t)
+    V = np.stack(v_t)
+    coords = {'lon':(['lon'],lon_c),'lat':(['lat'],lat_c),'time':(['time'],tempo)}    
+    ds = xr.Dataset({'u':(['time','lon','lat'],U),'v':(['time','lon','lat'],V)},coords=coords)
+   
+    ds['u']=ds.u.where(ds.u != 0,-9.8999995E15)
+    ds['v']=ds.v.where(ds.v != 0,-9.8999995E15)
+    
+    lon_attributtes = {'long_name': 'longitude',
+    'standard_name': 'longitude',\
+    'units': 'degrees_east',\
+    'valid_min': -180.0,\
+    'valid_max': 180.0}
+    
+    lat_attributtes = {'long_name': 'latitude',
+    'standard_name': 'latitude',\
+    'units': 'degrees_north',\
+    'valid_min': -90.0,\
+    'valid_max': 90.0}
+    time_atts = {'long_name':'time'}
+    ds.time.encoding['units'] = 'seconds since 1950-01-01 00:00:00'
+
+
+#    depth_attributtes = {'long_name': 'depth',
+#    'standard_name': 'depth',\
+#    'units': 'meters',\
+#    '_FillValue': -9.8999995E15,\
+#    'valid_min': np.min(self.grid['depth']),\
+#    'valid_max': np.max(self.grid['depth'])}
+   
+    ds.lon.attrs = lon_attributtes
+    ds.lat.attrs = lat_attributtes
+    #ds.depth.attrs = depth_attributtes
+    ds.time.attrs = time_atts
+   
+    encoding = {'u': {'_FillValue': -9.8999995E15},'v': {'_FillValue': -9.8999995E15}}
+    ds.to_netcdf(file_in[0:-5]+'.nc4',encoding=encoding)
+    ds.close()
+    return
+
+
+def download_HDF5_files(fecha):
+    # Fechas de inicio e fin
+    delta_time=timedelta(days=1)
+    data = fecha-delta_time
+    init_time=datetime(data.year,data.month,data.day)
+    end_time=datetime(data.year,data.month,data.day)
+    # Bucle que recorre las fechas.
+    date=init_time
+    print('-'*30)
+    while (date <= end_time):
+        print('Descargando archivo...')
+        print('-'*30) 
+        # Extraemos la info del objeto fecha en diferentes variables caracter
+        year=date.strftime("%Y")
+        month=date.strftime("%m")
+        day=date.strftime("%d")
+        # Definimos los nombres de los ficheros de entrada y salida
+        file_URL='https://mandeo.meteogalicia.es/thredds/fileServer/modelos/mohid/history/arousa/MOHID_Hydrodynamic_Arousa_%s_0000.hdf5' \
+        %(year+month+day)
+        file_name_out="MOHID_Arousa_%s_0000.hdf5" %(year+month+day)
+        try:
+            # Descargamos el fichero
+            with tqdm(unit='B', unit_scale=True, miniters=1, desc=file_name_out) as tqdm_instance:
+                urllib.request.urlretrieve(file_URL, file_name_out, reporthook=lambda block_num, block_size, total_size: 
+                tqdm_instance.update(block_num * block_size - tqdm_instance.n))
+        except urllib.error.HTTPError as e:
+            print(f'Error al descargar el archivo: {e}')
+            return None
+        else:
+            print(file_name_out)
+        finally:
+            # Avanzar la fecha
+            date=date+delta_time
+            print('-'*30)
+    return 'Done'
+
+if __name__ == '__main__':
+    fecha = datetime(2021,12,7,7)
+    download_HDF5_files(fecha)
+    ruta = os.getcwd()+'\\'
+    all_files = glob.glob(os.path.join(ruta +'*.hdf5'))
+    MOHID_hdf5_to_nc(all_files[0])
+    
+    
+    
+    
+    
+    
\ No newline at end of file
diff --git a/common/inout.py b/common/inout.py
new file mode 100644
index 0000000..ea0c72c
--- /dev/null
+++ b/common/inout.py
@@ -0,0 +1,23 @@
+import json
+from collections import OrderedDict
+
+
+def read_input(input_file, input_keys):
+    try:
+        with open(input_file, 'r') as f:
+            return json.load(f, object_pairs_hook=OrderedDict)
+    except FileNotFoundError:
+        print(f'File not found: {input_file} ')
+        if input('Do you want to create one (y/n)?') == 'n':
+            quit()
+
+        print(f'A {input_file} will be created with the next keys:\n')
+        json_obj = {}
+        for input_key in input_keys:
+            json_obj[input_key] = input(f'key: {input_key}?\n')
+        print(f'Writing a json_file: {input_file} with the next content:')
+        print(json.dumps(json_obj, indent=4))
+        with open(input_file, 'w') as json_file:
+            json.dump(json_obj, json_file, indent=4)
+        print('Done!\n')
+        return read_input(input_file, input_keys)
diff --git a/common/map.py b/common/map.py
new file mode 100644
index 0000000..c2d39c0
--- /dev/null
+++ b/common/map.py
@@ -0,0 +1,79 @@
+import cartopy.crs as ccrs
+import matplotlib.pyplot as plt
+
+
+class Map:
+
+    def __init__(self, projection=ccrs.PlateCarree()) -> None:
+        self.wms = None
+        self.layers = None
+        self.projection = projection
+        self.extension = None
+        self.ax = plt.axes(projection=self.projection)
+
+    def add_figure(self, figure_size=(14, 8)):
+        fig = plt.figure(figsize=figure_size)
+        self.ax = fig.add_subplot(1, 1, 1, projection=self.projection)
+        return fig
+
+    def set_extension(self, lat_south: float, lat_north: float, lon_west: float, lon_east: float) -> None:
+        self.extension = [lon_west, lon_east, lat_south, lat_north]
+        self.ax.set_extent(self.extension)
+    
+    def add_wms(self, wms: str, layers) -> None:
+        self.wms = wms
+        self.layers = layers
+        if type(layers) is list:
+            self.ax.add_wms(wms, layers)
+        elif type(layers) is str:
+            self.ax.add_wmts(wms, layers)
+
+
+    @staticmethod
+    def show():
+        plt.show()
+    
+
+def main():
+
+    arousa = Map()   
+    arousa.set_extension(41.7, 43.4, -9.5, -7.7)
+   
+    #arousa.add_wms('http://www.ign.es/wms-inspire/pnoa-ma?',
+    #['fondo', 'OI.OrthoimageCoverage'])
+
+    #arousa.add_wms('https://ideihm.covam.es/ihm-inspire/wms-regionesmarinas?', ['SR.Shoreline'])
+
+
+
+    #arousa.add_wms('https://wms.mapama.gob.es/sig/Agua/RiosCompPfafs/wms.aspx?', ['HY.PhysicalWaters.Waterbodies'])
+
+    arousa.add_wms('https://services.arcgisonline.com/arcgis/rest/services/Canvas/World_Light_Gray_Base/MapServer/WMTS?'
+                   , 'Canvas_World_Light_Gray_Base')
+    #arousa.add_wms('https://ows.emodnet-bathymetry.eu/wms?', ['emodnet:mean_atlas_land', 'coastlines'])
+    arousa.add_wms('https://ows.emodnet-bathymetry.eu/wms?', ['coastlines'])
+    #arousa.add_wms('https://map1c.vis.earthdata.nasa.gov/wmts-geo/wmts.cgi', 'VIIRS_CityLights_2012')
+    #arousa.add_wms('https://services.arcgisonline.com/arcgis/rest/services/Canvas/World_Dark_Gray_Base/MapServer/WMTS?',
+    #'Canvas_World_Dark_Gray_Base')
+    #arousa.add_wms("https://wms.mapama.gob.es/sig/EvaluacionAmbiental/Residuos/ParcelasAgric/wms.aspx?",
+                   #["Parcelas agrícolas de aplicación de lodos de EDAR"])
+    
+    arousa.add_wms("https://wms.mapama.gob.es/sig/Agua/EDAR/2021/wms.aspx?",
+                   ["Depuradoras de aguas residuales (Q2021. Dir 91/271/CEE)"])
+    arousa.add_wms("https://wms.mapama.gob.es/sig/Agua/PVERT/2021/wms.aspx?", 
+                   ["Puntos de vertido de depuradoras urbanas (Q2021. Dir 91/271/CEE)"])
+    arousa.add_wms("https://wms.mapama.gob.es/sig/Agua/CNV/wms.aspx?", 
+                   ["Censo Nacional de Vertidos"])
+    arousa.add_wms("https://wms.mapama.gob.es/sig/Agua/RedHidrograficaMDT/wms.aspx?", 
+                   ["Red hidrográfica básica MDT 100x100"])
+    
+    
+    #manager = plt.get_current_fig_manager()
+    #manager.window.showMaximized()
+    arousa.show()
+
+
+
+if __name__ == '__main__':
+    main()
+
diff --git a/common/readers/__init__.py b/common/readers/__init__.py
new file mode 100644
index 0000000..e69de29
diff --git a/common/readers/__pycache__/__init__.cpython-39.pyc b/common/readers/__pycache__/__init__.cpython-39.pyc
new file mode 100644
index 0000000..ee46a56
Binary files /dev/null and b/common/readers/__pycache__/__init__.cpython-39.pyc differ
diff --git a/common/readers/__pycache__/reader.cpython-39.pyc b/common/readers/__pycache__/reader.cpython-39.pyc
new file mode 100644
index 0000000..7fcacc9
Binary files /dev/null and b/common/readers/__pycache__/reader.cpython-39.pyc differ
diff --git a/common/readers/__pycache__/reader_HDF.cpython-39.pyc b/common/readers/__pycache__/reader_HDF.cpython-39.pyc
new file mode 100644
index 0000000..bbfe34e
Binary files /dev/null and b/common/readers/__pycache__/reader_HDF.cpython-39.pyc differ
diff --git a/common/readers/__pycache__/reader_NetCDF.cpython-39.pyc b/common/readers/__pycache__/reader_NetCDF.cpython-39.pyc
new file mode 100644
index 0000000..71c1722
Binary files /dev/null and b/common/readers/__pycache__/reader_NetCDF.cpython-39.pyc differ
diff --git a/common/readers/__pycache__/reader_factory.cpython-39.pyc b/common/readers/__pycache__/reader_factory.cpython-39.pyc
new file mode 100644
index 0000000..5236e5e
Binary files /dev/null and b/common/readers/__pycache__/reader_factory.cpython-39.pyc differ
diff --git a/common/readers/reader.py b/common/readers/reader.py
new file mode 100644
index 0000000..1baae7c
--- /dev/null
+++ b/common/readers/reader.py
@@ -0,0 +1,72 @@
+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+
+"""
+**reader.py**
+
+* *Purpose:* Abstract class of reader object. Reader object deals with reading HDF, NetCDF son classes
+
+* *python version:* 3.9
+* *author:* Pedro Montero
+* *license:* INTECMAR
+* *requires:*
+* *date:* 2021/10/13
+* *version:* 1.0.0
+* *date version* 2021/10/13
+
+
+"""
+
+from abc import ABC, abstractmethod
+
+
+class Reader(ABC):
+    """ Abstract class of Reader object"""
+    def __init__(self, file):
+        """ Open a file and get longitudes and latitudes and coordinates rank"""
+        self.file = file
+        self.dataset = self.open()
+        self.n_longitudes = None
+        self.n_latitudes = None
+        self.longitudes = self.get_longitudes()
+        self.latitudes = self.get_latitudes()
+        self.coordinates_rank = self.get_rank(self.longitudes)
+        self.ini_ntime = self.get_ini_ntime()
+        self.close()
+
+    @abstractmethod
+    def open(self):
+        pass
+
+    @abstractmethod
+    def close(self):
+        pass
+
+    @abstractmethod
+    def get_latitudes(self):
+        pass
+
+    @abstractmethod
+    def get_longitudes(self):
+        pass
+
+    @abstractmethod
+    def get_dates(self):
+        pass
+
+    @abstractmethod
+    def get_date(self, n_time):
+        pass
+
+    @abstractmethod
+    def get_ini_ntime(self):
+        pass
+
+    @staticmethod
+    def get_rank(array) -> int:
+        """ Number of dimensions of an array"""
+        return len(array.shape)
+
+
+
+
diff --git a/common/readers/reader_HDF.py b/common/readers/reader_HDF.py
new file mode 100644
index 0000000..ef23993
--- /dev/null
+++ b/common/readers/reader_HDF.py
@@ -0,0 +1,63 @@
+
+from datetime import datetime
+import numpy as np
+import h5py
+from .reader import Reader
+
+
+class ReaderHDF(Reader):
+
+    names = {
+        'northward_velocity': '/Results/velocity V/velocity V_',
+        'eastward_velocity': '/Results/velocity U/velocity U_'
+    }
+
+    def open(self):
+        self.dataset = h5py.File(self.file)
+        return self.dataset
+
+    def close(self):
+        self.dataset.close()
+
+    def get_latitudes(self):
+        lat_in = self.dataset['/Grid/Latitude']
+        if len(lat_in.shape) == 1:
+            self.n_latitudes = lat_in.shape[0]
+            return lat_in
+        elif len(lat_in.shape) == 2:
+            self.n_latitudes = lat_in.shape[1]
+            return lat_in[0, ]
+
+    def get_longitudes(self):
+        lon_in = self.dataset['/Grid/Longitude']
+        if len(lon_in.shape) == 1:
+            self.n_longitudes = lon_in.shape[0]
+            return lon_in
+        elif len(lon_in.shape) == 2:
+            self.n_longitudes = lon_in.shape[0]
+            return lon_in[:, 1]
+
+    def get_dates(self):
+        return self.dataset['/Time']
+
+    def get_date(self, n_time):
+
+        date_in = self.dataset['/Time/Time_' + str(n_time).zfill(5)]
+        return datetime(year=int(date_in[0]), month=int(date_in[1]), day=int(date_in[2]),
+                        hour=int(date_in[3]), minute=int(date_in[4]), second=int(date_in[5]))
+
+    def get_variable(self, name_var, n_time):
+        path = self.names[name_var]
+        variable = self.dataset[path + str(n_time).zfill(5)]
+        if len(variable.shape) == 2:
+            variable = np.transpose(variable)
+        elif len(variable.shape) == 3:
+            variable = np.transpose(variable, (0, 2, 1))
+        return variable
+
+    def get_ini_ntime(self):
+        return 1
+
+
+
+
diff --git a/common/readers/reader_NetCDF.py b/common/readers/reader_NetCDF.py
new file mode 100644
index 0000000..611f944
--- /dev/null
+++ b/common/readers/reader_NetCDF.py
@@ -0,0 +1,59 @@
+
+from datetime import datetime
+import netCDF4
+from .reader import Reader
+
+
+class ReaderNetCDF(Reader):
+
+    def open(self):
+        dataset = netCDF4.Dataset(self.file)
+        self.variables = dataset.variables
+        return dataset
+
+    def close(self):
+        self.dataset.close()
+
+    def get_latitudes(self):
+        lat_in = self.get_var('latitude')
+        if len(lat_in.shape) == 1:
+            self.n_latitudes = lat_in.shape[0]
+        elif len(lat_in.shape) == 2:
+            self.n_latitudes = lat_in.shape[1]
+        return lat_in
+
+    def get_longitudes(self):
+        lon_in = self.get_var('longitude')
+        if len(lon_in.shape) == 1:
+            self.n_longitudes = lon_in.shape[0]
+        elif len(lon_in.shape) == 2:
+            self.n_longitudes = lon_in.shape[0]
+        return lon_in
+
+    def get_dates(self):
+        times_in = self.get_var('time')
+        return netCDF4.num2date(times_in[:], units=times_in.units)
+
+    def get_date(self, n_time):
+        return self.get_dates()[n_time]
+
+    def get_variable(self, var_name, n_time):
+        return self.get_var(var_name)[n_time, ]
+
+    def get_var(self, var_name):
+        """Return values using the CF standard name, long_name or the simple name of a variable in a netCDF file."""
+        for var in self.variables:
+            for atributo in (self.variables[var].ncattrs()):
+                if atributo == 'standard_name':
+                    nome_atributo = (getattr(self.variables[var], 'standard_name'))
+                    if nome_atributo == var_name:
+                        return self.variables[var]
+                elif atributo == 'long_name':
+                    nome_atributo = (getattr(self.variables[var], 'long_name'))
+                    if nome_atributo == var_name:
+                        return self.variables[var]
+                elif var == var_name:
+                    return self.variables[var]
+
+    def get_ini_ntime(self):
+        return 0
diff --git a/common/readers/reader_factory.py b/common/readers/reader_factory.py
new file mode 100644
index 0000000..5738b27
--- /dev/null
+++ b/common/readers/reader_factory.py
@@ -0,0 +1,52 @@
+import sys
+from abc import ABC, abstractmethod
+from .reader import Reader
+from .reader_HDF import ReaderHDF
+from .reader_NetCDF import ReaderNetCDF
+
+
+class ReaderFactory(ABC):
+    """Basic class of factory readers"""
+
+    def __init__(self, file_in):
+        self.file_in = file_in
+
+    @abstractmethod
+    def get_reader(self) -> Reader:
+        """return a reader class"""
+
+
+class ReaderHDFFactory(ReaderFactory):
+    """Factory for Reader of HDF files"""
+
+    def get_reader(self) -> Reader:
+        return ReaderHDF(self.file_in)
+
+
+class ReaderNetCDFFactory(ReaderFactory):
+    """ Factoyr for Reader of NetCDF files"""
+
+    def get_reader(self) -> Reader:
+        return ReaderNetCDF(self.file_in)
+
+
+def read_factory(file_in) -> ReaderFactory:
+    """ Construct a reader factory based on the extension of file"""
+
+    factories = {
+        "nc": ReaderNetCDFFactory,
+        "nc4": ReaderNetCDFFactory,
+        "hdf": ReaderHDFFactory,
+        "hdf5": ReaderHDFFactory
+    }
+
+    extension = file_in.split('.')[-1]
+    if extension in factories:
+        return factories[extension](file_in)
+    print(f'Unknown file extension {extension}')
+    sys.exit(1)
+
+
+
+
+
diff --git a/common/test_readers.py b/common/test_readers.py
new file mode 100644
index 0000000..ee43b94
--- /dev/null
+++ b/common/test_readers.py
@@ -0,0 +1,27 @@
+
+from readers.reader_factory import read_factory
+
+
+def main():
+
+    factory = read_factory('../datos/drawmap/MOHID_Hydrodynamic_Vigo_20180711_0000.hdf5')
+    reader = factory.get_reader()
+    print(reader.longitudes)
+    print(reader.latitudes)
+    print(reader.get_date(1))
+    print(reader.get_date(10))
+
+    factory = read_factory('../datos/drawmap/wrf_arw_det1km_history_d05_20180711_0000.nc4')
+    reader = factory.get_reader()
+
+    print(reader.latitudes)
+    print(reader.longitudes)
+    print(reader.get_date(1))
+    print(reader.get_date(10))
+    print(f'number of longitudes = {reader.n_longitudes} and latitudes = {reader.n_latitudes}')
+    print(reader.get_variable('wind_module', 1))
+    reader.close()
+
+
+if __name__ == "__main__":
+    main()
diff --git a/common/vector.py b/common/vector.py
new file mode 100644
index 0000000..c47b0ec
--- /dev/null
+++ b/common/vector.py
@@ -0,0 +1,74 @@
+from math import pi, atan2, sin, cos
+import numpy as np
+
+
+class Vector:
+
+    def __init__(self, u=None, v=None, mod=None, theta=None, wind=False):
+        self.wind = wind
+
+        self.u = u
+        self.v = v
+        self.mod = mod
+        self.theta = theta
+
+        self.x = None
+        self.y = None
+
+        if self.mod is None or self.theta is None:
+            self.mod, self.theta = self.get_modtheta()
+
+        if self.u is None or self.v is None:
+            self.u, self.v = self.get_uv()
+
+    def set_point(self, x, y):
+        self.x = x
+        self.y = y
+
+    def get_modtheta(self):
+        """
+        From u,v velocity components return module and bearing of current
+
+        """
+        if self.wind:
+            coef = -1
+        else:
+            coef = 1
+        deg2rad = 180. / pi
+
+        mod = pow((pow(self.u, 2) + pow(self.v, 2)), .5)
+        theta = atan2(coef * self.u, coef * self.v) * deg2rad
+        return mod, theta
+
+    def get_uv(self):
+        """
+        From module, bearing of a current return u,v velocity components
+
+        """
+        if self.wind:
+            coef = -1
+        else:
+            coef = 1
+
+        rad2deg = pi / 180.
+        u = coef * self.mod * sin(self.theta * rad2deg)
+        v = coef * self.mod * cos(self.theta * rad2deg)
+        return u, v
+
+
+class NumpyVector(Vector):
+
+    def get_modtheta(self):
+        """
+        From u,v velocity components return module and bearing of current
+
+        """
+        if self.wind:
+            coef = -1
+        else:
+            coef = 1
+        deg2rad = 180. / pi
+
+        mod = pow((pow(self.u, 2) + pow(self.v, 2)), .5)
+        theta = np.arctan2(coef * self.u, coef * self.v) * deg2rad
+        return mod, theta
diff --git a/time90.py b/time90.py
new file mode 100644
index 0000000..3457645
--- /dev/null
+++ b/time90.py
@@ -0,0 +1,1339 @@
+# -*- coding: utf-8 -*-
+"""
+Created on Wed Mar  8 11:18:57 2023
+
+@author: josmo
+"""
+
+# Data processing
+import numpy as np
+# For additional mathematical functions.
+import math as m
+# For handling data in various formats
+import pandas as pd
+# For graphical representations
+import matplotlib.pyplot as plt
+# For dates and times
+from datetime import datetime, timedelta
+# For directories and system files
+import os
+# For secure file downloading
+import urllib
+# To show the progress of the download
+from tqdm import tqdm
+# For searching for files in a specific directory pattern
+import glob
+# For handling labeled multidimensional data
+import xarray as xr
+# For handling regular expressions
+import re
+# For handling images.
+from PIL import Image, ImageTk
+# To create graphical user interfaces (GUI).
+import tkinter as tk
+# To color printed text in the console
+from termcolor import colored
+# We import functions from other programs
+from common.readers.reader_factory import read_factory
+from common import read_input
+from common.boundarybox import BoundaryBox
+from common.drawcurrents import drawcurrents
+from common.map import Map
+from common.drawmap import*
+
+
+def distancia_esferica(lat1, lon1, lat2, lon2):
+    """
+    Calculates the distance between two points on a sphere given their latitude and longitude coordinates.
+    
+    Parameters:
+    -----------
+    lat1 : float
+        Latitude of the first point in degrees.
+    lon1 : float
+        Longitude of the first point in degrees.
+    lat2 : float
+        Latitude of the second point in degrees.
+    lon2 : float
+        Longitude of the second point in degrees.
+    
+    Returns:
+    --------
+    d : float
+        Distance between the two points in kilometers.
+    """
+    
+    r_tierra = 6371 # Radius of the Earth in kilometers
+    d_lat = np.radians(lat2 - lat1)
+    d_lon = np.radians(lon2 - lon1)
+    a = np.sin(d_lat/2)**2 + np.cos(np.radians(lat1))*np.cos(np.radians(lat2))*np.sin(d_lon/2)**2
+    #c = 2 * np.arctan2(np.sqrt(a), np.sqrt(1-a))
+    c = 2*np.arcsin(np.sqrt(a))
+    d = r_tierra * c # distance
+    return d
+
+
+def distancia_pitagoras(lat_input, lon_input, lats, lons):
+    """
+    Calculate the distance between two points using the Pythagorean theorem.
+    
+    Parameters:
+    -----------
+    lat_input: float
+        Latitude of the input point in degrees.
+    lon_input: float
+        Longitude of the input point in degrees.
+    lats: numpy.ndarray
+        Array of latitudes of the points to calculate the distance to, in degrees.
+    lons: numpy.ndarray
+        Array of longitudes of the points to calculate the distance to, in degrees.
+    
+    Returns:
+    --------
+    d: numpy.ndarray
+        Array of distances in kilometers between the input point and each of the points specified by `lats` and `lons`.
+
+    """
+    
+    R_t = 6371 # Radius of the Earth in kilometers
+    dy = (np.pi/180)*(lat_input-lats)*R_t
+    dx = (np.pi/180)*(lon_input-lons)*R_t*np.cos(np.radians(lats))
+    d = np.sqrt(dx**2+dy**2) # distance
+    #d = (2*np.pi*6371/360)*np.sqrt((lats-lat_input)**2 + (lons-lon_input)**2)
+    return d
+
+
+def distancia_min(lat_input, lon_input, lats, lons):
+    """
+    Calculates the minimum distance between an input point and a set of points defined by latitudes and longitudes.
+
+    Parameters
+    ----------
+    lat_input : float
+        Latitude of the input point in degrees.
+    lon_input : float
+        Longitude of the input point in degrees.
+    lats : array-like
+        Array of latitudes for the set of points in degrees.
+    lons : array-like
+        Array of longitudes for the set of points in degrees.
+
+    Returns
+    -------
+    index_min: int
+        Index of the station closest to the reference point
+    dist_min : float
+        Approximate distance in kilometers between the reference point and the closest station.
+        
+    If the function does not find a valid minimum distance (for example, if the reference point is too far away),
+    the distance and index obtained using the spherical distance method are returned.
+    """
+    
+    dist_esf = distancia_esferica(lat_input, lon_input, lats, lons)
+    dist_pit = distancia_pitagoras(lat_input, lon_input, lats, lons)
+    # Find the minimum distance and the corresponding index to identify the closest point
+    indice_min_esf = np.argmin(dist_esf)
+    indice_min_pit = np.argmin(dist_pit)
+    distancia_min_esf = dist_esf[indice_min_esf]
+    distancia_min_pit = dist_pit[indice_min_pit]
+    # If the minimum distances calculated by both methods are approximately the same and the indices of the points with the minimum distances are the same,
+    #then the function returns the index and the average of the two distances
+    if (indice_min_esf==indice_min_pit) and np.round(distancia_min_esf,0)==np.round(distancia_min_pit,0):
+        index_min = indice_min_esf
+        dist_min = np.round((distancia_min_esf+distancia_min_pit)/2,2)
+        return (index_min,dist_min)
+    else:
+        print('The spill is too far away from the observation stations.')
+        return (indice_min_esf,distancia_min_esf) # returns the index and distance calculated using the spherical distance method
+    
+
+def estaciones_proximas_rad_solar(Lat,Lon):
+    """
+    Given a latitude and longitude coordinate, this function finds the nearest weather observation station for 
+    solar radiation data, based on a pre-defined dataset of station locations. 
+   
+    Parameters
+    ----------
+    Lat : float
+        Latitude of the point of interest.
+    Lon : float
+        Longitude of the point of interest.
+
+    Returns
+    -------
+    index_rad: int
+        Index of the closest weather station in the MeteoGz_Rad dataframe.
+    dist_rad : float
+        distance between the point of interest and the closest weather observation
+        station, in kilometers.
+
+    """
+    
+    # Path to the data file
+    ruta = os.getcwd()+'\\common\\'
+    datafile = 'Rias_UWWTP.xlsx'
+    
+    # Read in the data
+    MeteoGz_Rad = pd.read_excel(ruta + datafile, sheet_name= 'MeteoGz_Rad')
+    Data_Links = pd.read_excel(ruta + datafile, sheet_name= 'Data_Links')
+    
+    # Coordinates of the MeteoGz weather stations for solar radiation data
+    MeteoGz_Rad_lats = MeteoGz_Rad.Lat.values
+    MeteoGz_Rad_lons = MeteoGz_Rad.Lon.values
+    
+    # Find the nearest station and its distance
+    index_rad,dist_rad = distancia_min(Lat, Lon, MeteoGz_Rad_lats, MeteoGz_Rad_lons)
+    
+    # Display relevant information
+    print('You can check the following links to obtain solar radiation data:\n')
+    print('- Solar Radiation Data from MeteoGalicia:\n', Data_Links.values[2][0])
+    print("\nThe nearest weather observation station to ({}, {}) ".format(Lat,Lon)+
+          "is {} ".format(MeteoGz_Rad.values[index_rad][0])+
+          "located at ({}, {}) and approximately {} km away.".format(
+              MeteoGz_Rad_lats[index_rad],MeteoGz_Rad_lons[index_rad],dist_rad))
+    print('-'*40)
+    return (index_rad, dist_rad)
+
+
+def estaciones_proximas_TS(Lat,Lon):
+    """
+    Find the closest buoy station to a given latitude and longitude, and return its index and distance.
+
+    Parameters
+    ----------
+    Lat : float
+        Latitude of the point of interest.
+    Lon : float
+        Longitude of the point of interest.
+
+    Returns
+    -------
+    index_TS: int
+        Index of the closest buoy station in the MeteoGZ_Buoys dataframe.
+    dist_TS : float
+        Distance (in km) between the input coordinates and the closest buoy station.
+
+    """
+    
+    # Path of the file of interest
+    path = os.getcwd()+'\\common\\'
+    datafile = 'Rias_UWWTP.xlsx'
+    
+    # Read the file
+    MeteoGZ_Buoys = pd.read_excel(path + datafile, sheet_name= 'MeteoGZ_Buoys')
+    Data_Links = pd.read_excel(path + datafile, sheet_name= 'Data_Links')
+
+    # Coordinates of the MeteoGz stations where to consult salt+temp
+    MeteoGZ_Buoys_lats = MeteoGZ_Buoys.Lat.values
+    MeteoGZ_Buoys_lons = MeteoGZ_Buoys.Lon.values
+    
+    # Get the closest position and distance
+    index_TS,dist_TS = distancia_min(Lat, Lon, MeteoGZ_Buoys_lats, MeteoGZ_Buoys_lons)
+    
+    # Show the relevant information
+    print('You can check the following links to obtain the T and S data:\n')
+    print('- T & S data from MeteoGalicia:\n', Data_Links.values[3][0])
+    # Show the result
+    print("\nThe nearest observation station to ({}, {}) ".format(Lat,Lon)+
+          "is {} ".format(MeteoGZ_Buoys.values[index_TS][0])+
+          "located at ({}, {}) at a distance of approximately {} km.".format(
+              MeteoGZ_Buoys_lats[index_TS],MeteoGZ_Buoys_lons[index_TS],dist_TS))
+    print('-'*40)
+    
+    print('\nYou can also check data from other links such as:')
+    print('- T & S data from the CTD profiles of INTECMAR:', Data_Links.values[4][0]) 
+    return(index_TS,dist_TS)
+
+
+def tryint(s):
+    """
+    Attempts to convert the input argument `s` into an integer. If it is not possible, returns `s` unchanged.
+    """
+    
+    try:
+        return int(s)
+    except ValueError:
+        return s
+
+
+def alphanum_key(s):
+    """
+    Turns a string `s` into a list of string and number chunks by splitting it at every number. Returns the resulting list.
+    
+    Example:
+    >>> alphanum_key("z23a")
+    ["z", 23, "a"]
+    """
+    
+    return [tryint(c) for c in re.split('([0-9]+)', s)]
+
+
+def human_sort(l):
+    """
+    Sorts a list `l` in the way that humans expect, by splitting each string element into chunks of strings and numbers
+    and sorting them according to the values of the chunks. The sorting is done in place.
+    """
+    
+    l.sort(key=alphanum_key)
+
+    
+def download_MOHID_files(fecha,key):
+    """
+    Downloads MOHID files for a given date and key
+
+    Parameters
+    ----------
+    fecha : datetime.datetime
+        Date for which to download files.
+    key : list
+        List with two elements, where the first element is the MOHID key and the second element is the file name.
+
+    Returns
+    -------
+    str
+        'Done' if all files were downloaded successfully, or None if there was an HTTP error during file download.
+        
+    """
+    
+    # Set start and end dates
+    delta_time = timedelta(days=1)
+    data = fecha - delta_time
+    init_time = datetime(data.year, data.month, data.day)
+    end_time = datetime(data.year, data.month, data.day)
+    
+    # Loop through the dates
+    date = init_time
+    print('-'*30)
+    while (date <= end_time):
+        print('Downloading file...')
+        print('-'*30) 
+        
+        # Extract year, month, and day from date object
+        year = date.strftime("%Y")
+        month = date.strftime("%m")
+        day = date.strftime("%d")
+        
+        # Define the file URLs and names
+        file_URL = 'http://mandeo.meteogalicia.gal/thredds/fileServer/mohid_' + key[0] + '/fmrc/files/%s/MOHID_' + key[1] + '_%s_0000.nc4' \
+                   %(year+month+day,year+month+day)
+        file_name_out = "MOHID_" + key[1] + "_%s_0000.nc4" %(year+month+day)
+
+        print(file_URL, file_name_out)
+        
+        try:
+            # Download the file
+            #urllib.request.urlretrieve(file_URL,file_name_out)
+            with tqdm(unit='B', unit_scale=True, miniters=1, desc=file_name_out) as tqdm_instance:
+                urllib.request.urlretrieve(file_URL, file_name_out, reporthook=lambda block_num, block_size, total_size: 
+                tqdm_instance.update(block_num * block_size - tqdm_instance.n))            
+        except urllib.error.HTTPError as e:
+            # If there is an HTTP error, print the error and return None
+            print(f' 02 Error downloading file: {e}')
+            return None
+        else:
+            # If the file was downloaded successfully, print a success message
+            print(file_name_out+' has been downloaded successfully.')
+        finally:
+            # Advance the date and print separator lines
+            date = date + delta_time
+            print('-'*30)
+    
+    # If all files were downloaded successfully, return 'Done'
+    return 'Done'
+
+
+def delete_nc4_files():
+    """
+    Deletes all .nc4 files in the current working directory.
+    This function searches for all files in the current working directory 
+    that have a .nc4 file extension and deletes them if they exist.
+
+    Returns:
+        None. 
+        The function does not return anything, it simply deletes the files.
+        
+    """
+    
+    # Get the current working directory path and create a 'ruta' variable with a backslash at the end
+    ruta = os.getcwd()+'\\'
+    
+    # Use the glob module to find all files with the ".nc4" extension in the 'ruta' directory
+    all_files = glob.glob(os.path.join(ruta +'*.nc4'))
+    
+    # Iterate over the list of files and delete them using os.remove()
+    for file_path in all_files:
+        if os.path.exists(file_path):
+            os.remove(file_path)
+            print(f"The file {file_path} has been successfully deleted.")
+        else:
+            print(f"The file {file_path} does not exist.")
+
+            
+def extract_TS_nc_file(time,lat,lon,depth):
+    """
+    Extracts temperature and salinity values from NetCDF files in the current working directory.
+
+
+    Parameters
+    ----------
+    time : str or datetime.datetime
+        Time at which temperature and salinity values are desired.
+    lat : float
+        Latitude at which temperature and salinity values are desired.
+    lon : float
+        Longitude at which temperature and salinity values are desired.
+    depth : float
+        Depth at which temperature and salinity values are desired.
+
+    Returns
+    -------
+    T: float
+        Temperature value.
+    S : float
+        Salinity value.
+
+    """
+    
+    # Get the current working directory and append a backslash
+    ruta = os.getcwd()+'\\'
+    # Find all files in the directory with a .nc4 extension
+    all_files = glob.glob(os.path.join(ruta +'*.nc4'))
+    # Sort the files in human-readable order
+    human_sort(all_files)
+
+    # Open the first file in the sorted list as an xarray dataset
+    ds = xr.open_dataset(all_files[0])
+    # Convert the time values in the dataset to pandas datetime objects
+    datas = pd.to_datetime(ds.time.values)
+
+    # Check if the provided time is in the dataset
+    if str(time) in datas.astype(str):
+        # If so, get the index of the corresponding time value
+        time_index = np.where(datas == time)[0][0]
+    else:
+        # If not, print an error message
+        print("The provided time is not in the dataset.")
+
+    # Find the index of the latitude value closest to the provided latitude
+    lat_index = np.abs(lat - ds.lat.values).argmin()
+    # Find the index of the longitude value closest to the provided longitude
+    lon_index = np.abs(lon - ds.lon.values).argmin()
+    # Find the index of the depth value closest to the provided depth
+    depth_index = np.abs(depth - ds.depth.values).argmin()
+
+    # Get the temperature value at the specified time, latitude, longitude, and depth
+    T = ds.temp.values[time_index][depth_index][lat_index][lon_index]
+    # Get the salinity value at the specified time, latitude, longitude, and depth
+    S = ds.salt.values[time_index][depth_index][lat_index][lon_index]
+
+    # Return a tuple of the temperature and salinity values
+    return (T, S)
+  
+def T90(T,S,iz):
+    """
+    Parameters
+    ----------
+    T : float
+        Surrounding temperature in ºC.
+    S : float
+        Surrounding water salinity (psu).
+    iz : float
+        Light radiation (W/m2) at deph z (m).
+
+    Returns
+    -------
+    T90 : float
+        T90 is the time in which 90% of E.Coli population is no longer detectable
+    """
+    
+    k = 2.533*(1.04**(T-20))*(1.012**S)+0.113*iz
+    T90 = (2.303/k)*24 #T90 in hours
+    return (T90)
+
+def start_T90():
+    """
+    Calculates the 90% mortality time of E.Coli by taking input from the user for Temperature (ºC), Salinity (psu), and 
+    Solar Radiation (W/m2) at a given depth z. The user needs to select the study Ria and provide the spill position 
+    (Lat, Lon), and the date in the format YYYY-MM-DD HH. 
+    
+    Returns None if the user selects to exit the program or if there is an error in the inputs.
+    
+    """
+    
+    # Prints the heading of the function
+    print(colored('\n - T90 CALCULATION: \n', attrs=['underline']))
+    #print('\n - T90 CALCULATION: \n')
+    
+    # Prompts the user for input and describes the function
+    print('This function aims to calculate the 90% mortality time of E.Coli.\n'
+      'For this, we will need the values of Temperature (ºC), Salinity (psu) and\n'
+      'Solar Radiation (W/m2) at a depth z.\n'
+      'The first step is to select the study Ria, provide the spill position (Lat, Lon), and the date:')
+    
+    # List of keywords to select the Ria
+    keyword = [['arousa','Arousa'],['vigo','Vigo'],['noiamuros','NoiaMuros','Noiamuros'],
+               ['artabro','Artabro']]
+    # Rias list
+    Rias = ['Ria de Arousa', 'Ria de Pontevedra-Vigo', 'Ria de Noia-Muros', 'Ria do Artabro']
+    
+    nam = None
+    while nam not in Rias:
+        # Displays the available options for Rias and prompts the user to select one
+        print('The available Rias for this study are:')
+        for i in Rias:
+            print(f'- {i}')
+        # Ask user for Ria selection
+        nam = input('Select one of the above options or type "exit" to exit: ')
+        if nam == 'exit':
+            return None
+            break
+        if nam not in Rias:
+            print(f'Error: "{nam}" is not a valid option. Please try again.')
+
+    # Select Ria based on user input
+    key = keyword[Rias.index(nam)]
+    # Initialize control variable
+    continuar = True
+    
+    # Prompts the user for latitude
+    while continuar:
+        try:
+            entrada = input("Enter the latitude in decimal format or type 'exit' to exit: ")
+            if entrada == 'exit':
+                continuar = False
+                break
+            else:
+                lat = float(entrada)
+                if -90 <= lat <= 90:
+                    #if not 42.383 <= lat <= 42.680:
+                        #print('La latitud no pertenece a la Ría de Arousa')
+                    break
+                else:
+                    print("Latitude must be between -90 and 90 degrees.")
+        except ValueError:
+            print("Please enter a valid decimal number.")
+    
+    # Prompts the user for longitude
+    while continuar:
+        try:
+            entrada = input("Enter the longitude in decimal format or type 'exit' to exit: ")
+            if entrada == 'exit':
+                continuar = False
+                break
+            else:
+                lon = float(entrada)
+                if -180 <= lon <= 180:
+                    break
+                else:
+                    print("Longitude must be between -180 and 180 degrees.")
+        except ValueError:
+            print("Please enter a valid decimal number.")
+    
+    # Ask user for date
+    while continuar:
+        try:
+            entrada = input("Enter the date in YYYY-MM-DD HH format or type 'exit' to exit: ")
+            if entrada == 'exit':
+                continuar = False
+                break
+            else:
+                fecha = datetime.strptime(entrada, '%Y-%m-%d %H')
+                break
+        except ValueError:
+            print("Please enter a valid date format (YYYY-MM-DD HH).")
+    
+    # Asks the user if they have the necessary data to calculate T90, or if they want to exit
+    while continuar:
+        try:
+            entrada = input("Do you know the data to calculate T90? [y/n] or type 'exit' to exit: ")
+            if entrada == 'exit':
+                continuar = False
+                break
+            elif entrada.lower() in ['y', 'n']:
+                answer1 = entrada.lower()
+                
+                # Ask for the value of z
+                while answer1=='y':
+                    z_input = input("Enter the value of the depth z (m) or type 'exit' to exit: ")
+                    if z_input == 'exit':
+                        answer1 = False
+                        continuar = False
+                        break
+                    else:
+                        try:
+                            z = float(z_input)
+                            break
+                        except ValueError:
+                            print("Please enter a valid number.")
+                            
+                # Ask for the value of T
+                while answer1=='y':
+                    T_input = input("Enter the value of temperature T (ºC) or type 'exit' to exit: ")
+                    if T_input == 'exit':
+                        answer1 = False
+                        continuar = False
+                        break
+                    else:
+                        try:
+                            T = float(T_input)
+                            break
+                        except ValueError:
+                            print("Please enter a valid number.")
+                            
+                # Ask for the value of S
+                while answer1=='y':
+                    S_input = input("Enter the value of salinity S (psu) or type 'exit' to exit: ")
+                    if S_input == 'exit':
+                        answer1 = False
+                        continuar = False
+                        break
+                    else:
+                        try:
+                            S = float(S_input)
+                            break
+                        except ValueError:
+                            print("Please enter a valid number.")
+                            
+                # Ask for the value of iz
+                while answer1=='y':
+                    iz_input = input("Enter the value of solar radiation iz (W/m2) or type 'exit' to exit: ")
+                    if iz_input == 'exit':
+                        answer1 = False
+                        continuar = False
+                        break
+                    else:
+                        try:
+                            iz = float(iz_input)
+                            break
+                        except ValueError:
+                            print("Please enter a valid number.")
+                            
+                while answer1=='n':
+                    try:
+                        entrada = input("Do you want the data to calculate T90 from observations[o] \nor numerical models[m]? [o/m] or type 'exit' to exit: ")
+                        if entrada == 'exit':
+                            answer1 = False
+                            continuar = False
+                            break
+                        elif entrada.lower() in ['o', 'm']:
+                            answer2 = entrada.lower()
+                            if answer2 == 'o':
+                                estaciones_proximas_rad_solar(lat,lon)
+                                estaciones_proximas_TS(lat,lon)
+                                
+                                # Ask for the value of z
+                                while answer2 == 'o':
+                                    z_input = input("Enter the value of depth z (m) or type 'exit' to exit: ")
+                                    if z_input == 'exit':
+                                        answer1 = False; answer2 = False
+                                        continuar = False
+                                        break
+                                    else:
+                                        try:
+                                            z = float(z_input)
+                                            break
+                                        except ValueError:
+                                            print("Please enter a valid number.")
+                                            
+                                # Ask for the value of T
+                                while answer2 == 'o':
+                                    T_input = input("Enter the value of temperature T (ºC) or type 'exit' to exit: ")
+                                    if T_input == 'exit':
+                                        answer1 = False; answer2 = False
+                                        continuar = False
+                                        break
+                                    else:
+                                        try:
+                                            T = float(T_input)
+                                            break
+                                        except ValueError:
+                                            print("Please enter a valid number.")
+                                            
+                                # Ask for the value of S
+                                while answer2 == 'o':
+                                    S_input = input("Enter the value of salinity S (psu) or type 'exit' to exit: ")
+                                    if S_input == 'exit':
+                                        answer1 = False; answer2 = False
+                                        continuar = False
+                                        break
+                                    else:
+                                        try:
+                                            S = float(S_input)
+                                            break
+                                        except ValueError:
+                                            print("Please enter a valid number.")
+                                            
+                                # Ask for the value of iz
+                                while answer2 == 'o':
+                                    iz_input = input("Enter the value of solar radiation iz (W/m2) or type 'exit' to exit: ")
+                                    if iz_input == 'exit':
+                                        answer1 = False; answer2 = False
+                                        continuar = False
+                                        break
+                                    else:
+                                        try:
+                                            iz = float(iz_input)
+                                            break
+                                        except ValueError:
+                                            print("Please enter a valid number.")
+                                            
+                            elif answer2 == 'm':
+                                print(colored("\nMeteoGalicia database:\n", attrs=['underline']))
+                                print("\nWe'll use MOHID model data from MeteoGalicia for T and S.")
+                                print("You can find the data at: http://mandeo.meteogalicia.gal/thredds/catalog.html")
+                                print("You'll have to look up the Solar Rad data.")
+                                estaciones_proximas_rad_solar(lat,lon)
+                                
+                                # Ask for the value of iz
+                                while answer2 == 'm':
+                                    iz_input = input("Enter the value of solar radiation iz (W/m2) or type 'exit' to exit: ")
+                                    if iz_input == 'exit':
+                                        answer1 = False; answer2 = False
+                                        continuar = False
+                                        break
+                                    else:
+                                        try:
+                                            iz = float(iz_input)
+                                            break
+                                        except ValueError:
+                                            print("Please enter a valid number.")
+                                            
+                                # Ask for the value of z
+                                while answer2 == 'm':
+                                    z_input = input("Enter the value of depth z (m) or type 'exit' to exit " +
+                                                    "(If you don't know it, enter 0 to calculate on the surface): ")
+                                    if z_input == 'exit':
+                                        answer1 = False; answer2 = False
+                                        continuar = False
+                                        break
+                                    else:
+                                        try:
+                                            z = float(z_input)
+                                            break
+                                        except ValueError:
+                                            print("Please enter a valid number.")
+                                            
+                                # Now, we download and perform automatic reading of the MOHID files.
+                                aux = download_MOHID_files(fecha,key) 
+                                if aux == 'Done':
+                                    # We extract the values of T and S
+                                    T,S = extract_TS_nc_file(fecha,lat,lon,z)
+                                    # We check if T or S are NaN
+                                    if m.isnan(T) or m.isnan(S):
+                                        print("Error: one of the variables is 'nan'.")
+                                        print("Try to enter values further offshore or shallower.")
+                                        answer1 = False; answer2 = False
+                                        continuar = False
+                                    else:
+                                        print("T and S values obtained correctly.")
+                                else:
+                                    break
+                                 
+                            break
+                        else:
+                            print("Please enter 'o' or 'm'.")
+                    except ValueError:
+                        print("Please enter a valid response.") 
+                break
+            else:
+                print("Please enter 'y' or 'n'.")
+        except ValueError:
+            print("Please enter a valid response.") 
+            
+
+    # Print user input values
+    if continuar:
+        print(colored("\nUser input values:\n", attrs=['underline']))
+        print("- Latitude:", lat)
+        print("- Longitude:", lon)
+        print("- Date:", fecha)
+        if answer1=='y' or answer2 == 'o' or aux == 'Done':
+            print("- T (ºC):", T)
+            print("- S (psu):", S)
+            print("- z (m):", z)
+            print("- iz (W/m2):", iz)
+            print(colored("\nCalculate T90:\n", attrs=['underline']))
+            # Calculate T90 and print the result
+            time_90 = T90(T,S,iz)
+            print("- The value of T90 is:", time_90, 'h')
+            
+    return(lat,lon,fecha,time_90,key)
+
+#-------------------LCS PART-----------------------------------
+def download_LCS_files(fecha,key):
+    """
+    This function downloads MYCOASTLCS netCDF files for a given date range and location.
+
+    Parameters
+    ----------
+    fecha : datetime
+        The date to start downloading files (datetime object).
+    key : list
+        The location key, containing two strings.
+
+    Returns
+    -------
+    str
+        A message indicating the function has finished ('Done') or None if there was an error.
+    """
+    
+    # Calculate start and end dates
+    delta_time = timedelta(days=1)
+    init_time = datetime(fecha.year, fecha.month, fecha.day) - timedelta(days=2)
+    end_time = datetime(fecha.year, fecha.month, fecha.day)
+
+    # Loop through each day in the date range
+    date = init_time
+    while (date <= end_time):
+        print('03:  Downloading file...')
+        print('-' * 30)
+        
+        # Extract date info into separate variables
+        year = date.strftime("%Y")
+        month = date.strftime("%m")
+        day = date.strftime("%d")
+        
+        # Remove location key if it is 'NoiaMuros'
+        if key[1] == 'NoiaMuros':
+            key.pop(1)
+            
+        # Define URL and file name based on date and location key
+        #http: // thredds - gfnl.usc.es / thredds / fileServer / MYCOASTLCS / MYCOASTLCS_Vigo_20230318.nc
+        file_URL = 'http://thredds-gfnl.usc.es/thredds/fileServer/MYCOASTLCS/MYCOASTLCS_' + key[1] + '_%s.nc' \
+                   % (year + month + day)
+        file_name_out = "MYCOASTLCS_" + key[1] + "_%s.nc" % (year + month + day)
+
+        print(file_URL, file_name_out)
+
+        try:
+            # Download the file
+            with tqdm(unit='B', unit_scale=True, miniters=1, desc=file_name_out) as tqdm_instance:
+                urllib.request.urlretrieve(file_URL, file_name_out, reporthook=lambda block_num, block_size, total_size:
+                tqdm_instance.update(block_num * block_size - tqdm_instance.n))
+        except urllib.error.HTTPError as e:
+            print(f'Error downloading file: {e}')
+            return None
+        else:
+            print(file_name_out + ' downloaded successfully.')
+        finally:
+            # Advance the date
+            date = date + delta_time
+            print('-' * 30)
+            
+    return 'Done'
+
+
+def delete_nc_files():
+    """
+    This function deletes all netCDF files in the current working directory.
+
+    Returns:
+        None
+    """
+    
+    # Get current working directory and find all netCDF files
+    ruta = os.getcwd()+'\\'
+    all_files = glob.glob(os.path.join(ruta +'*.nc'))
+
+    # Loop through files and delete them if they exist
+    for file_path in all_files:
+        if os.path.exists(file_path):
+            os.remove(file_path)
+            print(f"File {file_path} deleted successfully.")
+        else:
+            print(f"File {file_path} does not exist.")
+            
+
+def delete_png_files(ruta):
+    """
+    This function deletes all PNG files in the specified directory.
+
+    Parameters
+    ----------
+    ruta : str
+        The path to the directory containing the PNG files.
+
+    Returns
+    -------
+    None.
+
+    """
+    
+    # Find all PNG files in directory
+    all_files = glob.glob(os.path.join(ruta +'*.png'))
+
+    # Loop through files and delete them if they exist
+    for file_path in all_files:
+        if os.path.exists(file_path):
+            os.remove(file_path)
+            print(f"File {file_path} deleted successfully.")
+        else:
+            print(f"File {file_path} does not exist.")
+            
+def delete_gif_files():
+    """
+    This function deletes all GIF files in the current working directory.
+
+    Returns:
+        None
+    """
+    
+    # Get current working directory and find all GIF files
+    ruta = os.getcwd()+'\\'
+    all_files = glob.glob(os.path.join(ruta +'*.gif'))
+
+    # Loop through files and delete them if they exist
+    for file_path in all_files:
+        if os.path.exists(file_path):
+            os.remove(file_path)
+            print(f"File {file_path} deleted successfully.")
+        else:
+            print(f"File {file_path} does not exist.")
+            
+            
+def extract_time_LCS_file(time):
+    """
+    This function extracts the time index and file index of the nearest time in LCS netCDF files.
+
+    Parameters
+    ----------
+    time : str
+        A string representing the time to search for, in format YYYY-MM-DDTHH:MM:SS.
+
+    Returns
+    -------
+    list or None
+        A list containing the time index and file index if the time is found, or None if not found.
+
+    """
+    
+    # Get current working directory and find all netCDF files
+    ruta = os.getcwd()+'\\'
+    all_files = glob.glob(os.path.join(ruta +'*.nc'))
+    human_sort(all_files)
+
+    # Open the first two files and get their time values
+    ds1 = xr.open_dataset(all_files[0])
+    ds2 = xr.open_dataset(all_files[1])
+    datas1 = pd.to_datetime(ds1.time.values)
+    datas2 = pd.to_datetime(ds2.time.values)
+
+    # Check if the time is in either of the files
+    if str(time) in datas1.astype(str):
+        time_index = np.where(datas1 == time)[0][0]
+        file_index = 0
+        return [time_index,file_index]
+    elif str(time) in datas2.astype(str):
+        time_index = np.where(datas2 == time)[0][0]
+        file_index = 1
+        return [time_index,file_index]
+    else:
+        print("The date is not in the files.")
+        return None
+    
+
+def read_inputs(input_file):
+    """Read keywords for options"""
+    input_keys = ['path_in',
+                  'file_in',
+                  'path_out'
+                  'file_out',
+                  'nx',
+                  'ny',
+                  'resolution',
+                  'scale',
+                  'n_time',
+                  'n_level',
+                  'title',
+                  'style',
+                  'limits',
+                  'vector',
+                  'scalar',
+                  'wms_url',
+                  'wms_layers']
+    return read_input(input_file, input_keys)
+
+def crear_gift(ruta_maps, ruta_gift, vel_maps, loop):
+    """
+    Creates an animated GIF file from a directory of PNG images.
+
+    Parameters
+    ----------
+    ruta_maps : str
+        Path to the directory containing the PNG images.
+    ruta_gift : str
+        Path to the output GIF file.
+    vel_maps : float
+        Speed of the animation in frames per second.
+    loop : int
+        Number of loops for the animation (0 for infinite loop).
+
+    Returns
+    -------
+    None.
+
+    """
+    
+    # Set the path to the directory of images
+    image_directory = ruta_maps
+    
+    # Get the list of file names in the directory
+    image_files = os.listdir(image_directory)
+    
+    # Sort the list of file names alphabetically
+    image_files = sorted(image_files)
+    
+    # Create a list of Image objects from Pillow
+    image_list = []
+    for filename in image_files:
+        if filename.endswith(".png"):
+            image_path = os.path.join(image_directory, filename)
+            image = Image.open(image_path)
+            image_list.append(image)
+    
+    # Save the list of images as an animated GIF file
+    gif_path = ruta_gift
+    frame_duration = 1e3 / vel_maps # Duration of each frame in milliseconds
+    image_list[0].save(gif_path, save_all=True, append_images=image_list[1:],
+                       duration=frame_duration, loop=loop)
+
+
+def mostrar_gift(ruta_gift):
+    """
+    Displays an animated GIF file in a new window.
+
+    Parameters
+    ----------
+    ruta_gift : str
+        The path to the GIF file to be displayed.
+
+    Returns
+    -------
+    None.
+
+    """
+    
+    # Path of the GIF file to be displayed
+    gif_path = ruta_gift
+    
+    # Open the GIF file using Pillow
+    gif_im = Image.open(gif_path)
+    
+    # Create a new window using Tkinter
+    root = tk.Tk()
+    canvas = tk.Canvas(root, width=gif_im.width, height=gif_im.height)
+    canvas.pack()
+    
+    # Convert each frame of the animation into a PhotoImage object
+    frames = []
+    for frame in range(0, gif_im.n_frames):
+        gif_im.seek(frame)
+        frame_im = ImageTk.PhotoImage(gif_im)
+        frames.append(frame_im)
+    
+    
+    def update_frame(frame_number):
+        """
+        Update the canvas with each frame of the animation.
+
+        Parameters
+        ----------
+        frame_number : int
+            The current frame number.
+
+        Returns
+        -------
+        None.
+
+        """
+        # Delete all items on the canvas
+        canvas.delete("all")
+        # Add the current frame to the canvas
+        print(frame_number, frames[frame_number])
+        canvas.create_image(0, 0, image=frames[frame_number], anchor="nw")
+        # Call the update_frame function again after a delay of 500 milliseconds
+        root.after(500, update_frame, (frame_number + 1) % len(frames))
+    
+    # Start the animation
+    update_frame(0)
+    
+    # Run the tkinter loop
+    root.mainloop()
+    
+def main_drawmap(inputs,start,end,ref_lon,ref_lat,key,draw_scale):
+    """
+    Creates maps based on user inputs.
+
+    Parameters
+    ----------
+    inputs : dict
+        A dictionary containing user-defined inputs.
+    start : int
+        The start date index of the data to be plotted.
+    end : int
+        The end date index of the data to be plotted.
+    ref_lon : float
+        A reference longitude for the plot.
+    ref_lat : float
+        A reference latitude for the plot.
+    key : list
+        The location key, containing two strings.
+    draw_scale : str
+        Parameter to change the color scale of the map. If you enter 'jet' use the jet scale. 
+        Otherwise it will use the viridis scale by default.
+
+    Returns
+    -------
+    None.
+
+    """
+    
+    # Start
+    print(colored("\nCreating maps:\n", attrs=['underline']))
+    
+    if inputs['n_time'] == 'all':
+        # Draw map for start-end period
+        draw_map_24(inputs,start,end,ref_lon, ref_lat,key,draw_scale)
+    else:
+        if inputs['vector']:
+            # Draw vector map
+            draw_map_vector(inputs, inputs['n_time'])
+        elif inputs['scalar']:
+            # Draw scalar map
+            draw_map_scalar(inputs, inputs['n_time'])
+        else:
+            # Print error message and quit
+            print('You must choose between scalar and vector maps in the drawmap.json file.\nWords: vector, scalar')
+            quit(1)
+            
+def delete_all():
+    """
+    Deletes all generated files.
+
+    Returns
+    -------
+    None.
+
+    """
+    
+    # Deletes all NetCDF files
+    delete_nc_files()
+    
+    # Deletes all NetCDF4 files
+    delete_nc4_files()
+    
+    # Deletes all PNG files in the maps directory
+    delete_png_files(os.getcwd()+'\\maps\\')
+    
+    # Deletes all PNG files in the current directory
+    delete_png_files(os.getcwd()+'\\')
+    
+    # Deletes all GIF files
+    delete_gif_files()
+    
+    return None
+
+def jet_map(start_index,file_index,T9,nam,key,ref_lon, ref_lat):
+    """
+    Creates maps based on user inputs.
+
+    Parameters
+    ----------
+    start_index : int
+        The start date index of the data to be plotted.
+    file_index : int
+        The index of the start file of the data to be plotted.
+    T9 : int
+        Time duration to be plotted.
+    nam : str
+        Name of the variable to represent
+    key : list
+        The location key, containing two strings.
+    ref_lon : float
+        A reference longitude for the plot.
+    ref_lat : float
+        A reference latitude for the plot.
+
+    Returns
+    -------
+    None.
+
+    """
+    
+    # Get all netCDF files in the current directory
+    all_files = glob.glob(os.path.join('*.nc'))
+
+    # Set the path for the drawmap.json file
+    json_path = os.getcwd() + '\\common\\'
+
+    # Read the input values from the drawmap.json file
+    inputs = read_inputs(json_path + 'drawmap.json')
+    inputs['path_in'] = os.getcwd()
+    inputs['path_out'] = os.getcwd()
+    inputs['file_out'] = "MYCOAST_" + nam + "_" + key[1] + "_during_T90"
+    inputs['scalar_magnitude'] = nam
+    inputs['title'] = "MYCOAST-" + nam + " Ria " + key[1] + ":"
+
+    # Open the netCDF file for the specified time range
+    time_nc = xr.open_dataset(all_files[file_index]).isel(time=slice(start_index, (start_index+1)))
+
+    # Open the netCDF file for the entire day
+    ds_aux = xr.open_dataset(all_files[file_index])
+
+    # Calculate the sum of the variable values over the specified time range
+    if (start_index+T9)>=24:
+            end = 24
+            aux = (ds_aux[nam]).isel(time=slice(start_index, end)).sum(dim='time', skipna=False).values
+            T9 = T9-(end-start_index)
+            if T9 == 0:
+                time_nc[nam][0] = aux
+                time_nc.to_netcdf('condensed.nc')
+            elif  0<T9<=24:
+                ds_aux2 = xr.open_dataset(all_files[file_index+1])
+                aux2 = (ds_aux2[nam]).isel(time=slice(0, T9)).sum(dim='time', skipna=False).values
+                time_nc[nam][0] = aux+aux2
+                time_nc.to_netcdf('condensed.nc')
+            elif T9>24:
+                start_index = 0; end = 24
+                ds_aux2 = xr.open_dataset(all_files[file_index+1])
+                aux2 = (ds_aux2[nam]).isel(time=slice(start_index, end)).sum(dim='time', skipna=False).values
+                T9 = T9-(end-start_index)
+                ds_aux3 = xr.open_dataset(all_files[file_index+2])
+                aux3 = (ds_aux3[nam]).isel(time=slice(0, T9)).sum(dim='time', skipna=False).values
+                time_nc[nam][0] = aux+aux2+aux3
+                time_nc.to_netcdf('condensed.nc')
+    else:
+            end = start_index+T9
+            aux = (ds_aux[nam]).isel(time=slice(start_index, end)).sum(dim='time', skipna=False).values
+            time_nc[nam][0] = aux
+            time_nc.to_netcdf('condensed.nc')  
+    
+    # Set the input file name for main_drawmap function
+    inputs['file_in'] = 'condensed.nc'
+    # Set the scale for the map to 'jet'
+    draw_scale = 'jet'
+
+    # Create the jet map
+    main_drawmap(inputs,0,1,ref_lon, ref_lat,key,draw_scale)
+
+    # Directory where the images are located
+    directory = os.getcwd()
+    # Get the list of files in the directory
+    files = os.listdir(directory)
+    # Filter to get the name of the file that meets the conditions
+    file_name = next((file for file in files if file.startswith(inputs['file_out'])), None)
+    if file_name is not None:
+        # Full path of the file
+        file_path = os.path.join(directory, file_name)
+        # Open the image
+        image = Image.open(file_path)
+        # Show the image
+        image.show()
+    else:
+        print("No image found with the specified name.")
+
+
+def start_MYCOAST(fecha,T90,ref_lon, ref_lat,key):
+    """
+    Downloads MYCOAST data for a given date and key, and creates an animation
+    of the data for a given time window using the `main_drawmap` function.
+
+    Parameters
+    ----------
+    fecha : str
+        A string representing the date in the format YYYYMMDD.
+    T90 : float
+        A float representing the time window (in hours) for which to create the animation.
+    ref_lon : float
+        A reference longitude for the plot.
+    ref_lat : float
+        A reference latitude for the plot.
+    key : list
+        The location key, containing two strings.
+
+    Returns
+    -------
+    None.
+
+    """
+    
+    print(colored("\nMYCOAST_LCS Analysis:\n", attrs=['underline']))
+
+    # Download the MYCOAST data for the specified date and variable
+    aux = download_LCS_files(fecha,key)
+    if aux == 'Done':
+        # Get all netCDF files in the current directory
+        all_files = glob.glob(os.path.join('*.nc'))
+
+        # Extract the start index and file index for the specified time range
+        T9 = int(np.round(T90,0))
+        index = extract_time_LCS_file(fecha)
+        start = index[0]
+        file_index = index[1]
+
+        # Get the list of variables in the netCDF file
+        dt = xr.open_dataset(all_files[0])
+        long_names = []
+        for i in dt.data_vars:
+            long_names.append(i)
+        
+        nam = None
+        while nam not in long_names:
+            # Prompt the user to select the variable to plot
+            print('The available variables are:')
+            for ln in long_names:
+                print(f'- {ln}')
+            nam = input('Enter the long_name of the variable to plot or type "exit" to quit: ')
+            if nam == 'exit':
+                break
+            if nam not in long_names:
+                print(f'Error: "{nam}" is not a valid long_name. Please try again.')
+        
+        # Create a plot for the specified time range
+        jet_map(start,file_index,T9,nam,key,ref_lon, ref_lat)
+
+        # Set the input parameters for the plot animation
+        json_path = os.getcwd() + '\\common\\'
+        inputs = read_inputs(json_path + 'drawmap.json')
+        inputs['path_in'] = os.getcwd()
+        inputs['path_out'] = os.getcwd()+'\\maps'
+        inputs['file_out'] = "MYCOAST_" + nam + "_" + key[1]
+        inputs['scalar_magnitude'] = nam
+        inputs['title'] = "MYCOAST-" + nam + " Ria " + key[1] + ":"
+        
+        # Create a plot for each time window within the specified time range
+        for i in all_files:
+            inputs['file_in'] = i
+          
+            if (start+T9)>=24:
+                end = 24
+            else:
+                end = start+T9
+            
+            # Set the scale for the map to 'viridis' (default)
+            draw_scale = 'viridis'
+            main_drawmap(inputs,start,end,ref_lon, ref_lat,key,draw_scale)
+            
+            T9 = T9-(end-start)
+            if T9 <= 0:
+                break
+            else:
+                start = 0
+            
+        # Create a GIF animation of the plots
+        # Path where the maps are stored
+        ruta_maps = os.getcwd()+'\\maps'
+        # Path to store the generated GIF animation
+        ruta_gift = os.getcwd() + '\\animacion.gif'
+        # Create the gift with a duration of 0.5 seconds each frame and in an infinite loop
+        crear_gift(ruta_maps, ruta_gift, 2, 0)
+        # Display the generated GIF animation
+        mostrar_gift(ruta_gift)
+             
+if __name__ == '__main__':
+    # Delete any existing files from previous executions
+    delete_all()
+    # Start the T90 calculation and user input process
+    lat,lon,fecha,time_90,key = start_T90()
+    # Start the MYCOAST analysis with the user inputs from the previous step
+    start_MYCOAST(fecha,time_90,lon,lat,key)