Repository: meteocat/pymica (language: python)

Issue:
Remove web browser when creating clusters
One feature of the current code is that clusters obtained by a Kmeans algorithm can be changed by the user through a web page. This is a bit problematic and it would be better to just return the cluster files when they are created.
Remove web browser when creating clusters
One feature of the current code is that clusters obtained by a Kmeans algorithm can be changed by the user through a web page. This is a bit problematic and it would be better to just return the cluster files when they are created.
Remove web browser when creating clusters
One feature of the current code is that clusters obtained by a Kmeans algorithm can be changed by the user through a web page. This is a bit problematic and it would be better to just return the cluster files when they are created.


Discussion / hints:

Fix patch:
diff --git a/pymica/utils/clusters.py b/pymica/utils/clusters.py
new file mode 100644
index 0000000..ce7cfcc
--- /dev/null
+++ b/pymica/utils/clusters.py
@@ -0,0 +1,176 @@
+"""Tool for creating clusters with the K-Means method.
+"""
+import json
+from os.path import dirname
+
+import numpy as np
+from osgeo import gdal, ogr, osr
+from scipy.ndimage import gaussian_filter
+from sklearn.cluster import KMeans
+
+from pymica.utils.geotools import get_utm_epsg_from_lonlat, reproject_point
+
+
+def create_clusters(locations: list | str, n_clusters: int, geojson_file: str) -> None:
+    """Group locations in a specified number of clusters using K-Means algorithm.
+
+    Args:
+        locations (list | str): List of dictionaries including location data or the
+                                file path with data as json format. Data for each
+                                location is a dictionary with 'id', 'alt', 'lon' and
+                                'lat' as keys.
+        n_clusters (int): Number of clusters to create.
+        geojson_file (str): File path where clusters are saved.
+    """
+    if isinstance(locations, str):
+        f_p = open(locations)
+        locations = json.load(f_p)
+
+    positions = np.zeros([len(locations), 2])
+    utm_proj4 = get_utm_epsg_from_lonlat(locations[0]["lon"], locations[0]["lat"])
+
+    positions_list = []
+    for i, loc_value in enumerate(locations):
+        positions[i][0] = loc_value["lon"]
+        positions[i][1] = loc_value["lat"]
+
+        utm_coordinates = reproject_point(
+            (loc_value["lon"], loc_value["lat"]), "EPSG:4326", utm_proj4
+        )
+        positions_list.append(utm_coordinates)
+    positions_list = np.array(positions_list)
+
+    kmeans = KMeans(n_clusters=n_clusters, random_state=0).fit(positions_list)
+
+    out_geojson = {"type": "FeatureCollection", "features": []}
+    for i, labels_value in enumerate(kmeans.labels_):
+        locations[i]["cluster"] = labels_value
+        out_geojson["features"].append(
+            {
+                "type": "Feature",
+                "geometry": {
+                    "type": "Point",
+                    "coordinates": [locations[i]["lon"], locations[i]["lat"]],
+                },
+                "properties": {
+                    "cluster": int(labels_value),
+                    "id": locations[i]["id"],
+                    "alt": locations[i]["alt"],
+                },
+            }
+        )
+
+    try:
+        with open(geojson_file, "w") as f_p:
+            json.dump(out_geojson, f_p)
+    except FileNotFoundError:
+        raise FileNotFoundError(dirname(geojson_file) + " directory does not exist.")
+
+
+def create_reprojected_geometries(geojson_file: str, epsg: int) -> ogr.DataSource:
+    """Reprojection of an ogr file to specified projection.
+    Taken from: https://pcjericks.github.io/py-gdalogr-cookbook/projection.html
+
+    Args:
+        geojson_file (str): Path to a GeoJSON file with geometries.
+        epsg (int): EPSG code of the new projection
+
+    Raises:
+        ValueError: If EPSG code is wrong.
+        IOError: If `geojson_file` does not exist.
+
+    Returns:
+        ogr.DataSource: Reprojected geometries.
+    """
+    out_proj = osr.SpatialReference()
+    out_proj.ImportFromEPSG(epsg)
+    if len(out_proj.ExportToPrettyWkt()) <= 1:
+        raise ValueError("Wrong EPSG code: {}".format(epsg))
+
+    in_proj = osr.SpatialReference()
+    in_proj.ImportFromEPSG(4326)
+
+    transf = osr.CoordinateTransformation(in_proj, out_proj)
+
+    in_ds = ogr.Open(geojson_file)
+    if in_ds is None:
+        raise IOError("File {} doesn't exist".format(geojson_file))
+    in_layer = in_ds.GetLayer()
+
+    mem_driver = ogr.GetDriverByName("MEMORY")
+    proj_ds = mem_driver.CreateDataSource("memData")
+    proj_layer = proj_ds.CreateLayer(
+        "clusters", out_proj, geom_type=ogr.wkbMultiPolygon
+    )
+
+    in_layer_def = in_layer.GetLayerDefn()
+    for i in range(0, in_layer_def.GetFieldCount()):
+        field_def = in_layer_def.GetFieldDefn(i)
+        proj_layer.CreateField(field_def)
+
+    out_layer_def = proj_layer.GetLayerDefn()
+
+    feature = in_layer.GetNextFeature()
+    while feature:
+        geom = feature.GetGeometryRef()
+        geom.Transform(transf)
+        proj_feat = ogr.Feature(out_layer_def)
+        proj_feat.SetGeometry(geom)
+        for i in range(0, out_layer_def.GetFieldCount()):
+            proj_feat.SetField(
+                out_layer_def.GetFieldDefn(i).GetNameRef(), feature.GetField(i)
+            )
+        proj_layer.CreateFeature(proj_feat)
+        proj_feat = None
+        feature = in_layer.GetNextFeature()
+
+    return proj_ds
+
+
+def rasterize_clusters(
+    ds_in: ogr.DataSource, raster_config: dict, sigma: float = 15
+) -> None:
+    """Rasterize clusters from the GeoJSON file generated by
+
+    Args:
+        ds_in (ogr.DataSource): Input layer.
+        raster_config (dict): The output raster characteristics. Must include
+            'out_file' (path to save the rasterized clusters), 'size'
+            (raster size (x, y)) and 'geotransform' ([ul_x, x_res, x_rot,
+            ul_y, y_rot, y_res]) keys.
+        sigma (float, optional): Sigma parameter for a Gaussian filter. Defaults to 15.
+
+    Raises:
+        ValueError: If `raster_config` does not have all required keys.
+    """
+    if raster_config.keys() < {"out_file", "size", "geotransform"}:
+        raise KeyError(
+            "`out_file`, `size`, `geotransform` must be in the `raster_config` "
+            "parameter."
+        )
+
+    layer = ds_in.GetLayer()
+    num_layers = layer.GetFeatureCount()
+    proj = layer.GetSpatialRef()
+
+    driver = gdal.GetDriverByName("GTIFF")
+    ds_out = driver.Create(
+        raster_config["out_file"],
+        raster_config["size"][0],
+        raster_config["size"][1],
+        num_layers,
+        gdal.GDT_Float32,
+    )
+    ds_out.SetGeoTransform(raster_config["geotransform"])
+    ds_out.SetProjection(proj.ExportToWkt())
+
+    for i in range(num_layers):
+        layer.SetAttributeFilter("cluster={}".format(i))
+        gdal.RasterizeLayer(ds_out, [i + 1], layer, burn_values=[1])
+
+    data = ds_out.ReadAsArray().astype(np.float32)
+    for i in range(num_layers):
+        blurred = gaussian_filter(data[i], sigma)
+        ds_out.GetRasterBand(i + 1).WriteArray(blurred)
+
+    ds_out = None
diff --git a/pymica/utils/create_cluster_files.py b/pymica/utils/create_cluster_files.py
deleted file mode 100644
index fec31c5..0000000
--- a/pymica/utils/create_cluster_files.py
+++ /dev/null
@@ -1,102 +0,0 @@
-'''
-Generates the cluster rasters from the geometries.
-The generated file can be used to calculate the multiregressions
-using clusters.
-'''
-from osgeo import ogr, osr, gdal
-import numpy
-from scipy.ndimage import gaussian_filter
-
-
-def create_reprojected_geoms(file_name, epsg):
-    '''Reprojects an ogr file to the desired projection
-    Taken from:
-    https://pcjericks.github.io/py-gdalogr-cookbook/projection.html#reproject-a-layer
-
-    Args:
-        file_name (str): The file to reproject path
-        epsg (int): The new projection EPSG code
-
-    Returns:
-        osgeo.ogr.DataSource: The reprojected ogr datasource object
-    '''
-    out_proj = osr.SpatialReference()
-    out_proj.ImportFromEPSG(epsg)
-    if len(out_proj.ExportToPrettyWkt()) <= 1:
-        raise ValueError("Wrong EPSG code: {}".format(epsg))
-
-    in_proj = osr.SpatialReference()
-    in_proj.ImportFromEPSG(4326)
-
-    transf = osr.CoordinateTransformation(in_proj, out_proj)
-
-    in_ds = ogr.Open(file_name)
-    if in_ds is None:
-        raise IOError("File {} doesn't exist".format(file_name))
-    in_layer = in_ds.GetLayer()
-
-    mem_driver = ogr.GetDriverByName('MEMORY')
-    proj_ds = mem_driver.CreateDataSource('memData')
-    proj_layer = proj_ds.CreateLayer("clusters", out_proj,
-                                     geom_type=ogr.wkbMultiPolygon)
-
-    in_layer_def = in_layer.GetLayerDefn()
-    for i in range(0, in_layer_def.GetFieldCount()):
-        field_def = in_layer_def.GetFieldDefn(i)
-        proj_layer.CreateField(field_def)
-
-    out_layer_def = proj_layer.GetLayerDefn()
-
-    feature = in_layer.GetNextFeature()
-    while feature:
-        geom = feature.GetGeometryRef()
-        geom.Transform(transf)
-        proj_feat = ogr.Feature(out_layer_def)
-        proj_feat.SetGeometry(geom)
-        for i in range(0, out_layer_def.GetFieldCount()):
-            proj_feat.SetField(out_layer_def.GetFieldDefn(i).GetNameRef(),
-                               feature.GetField(i))
-        proj_layer.CreateFeature(proj_feat)
-        proj_feat = None
-        feature = in_layer.GetNextFeature()
-
-    return proj_ds
-
-
-def rasterize_clusters(ds_in, out_conf, sigma=15):
-    """Takes the clusters file the GeoJSON generated by the web interface
-    and rasterizes and blurs them so they can be used in the pymica functions.
-    The output projection is the same as the input layers
-
-    Args:
-        ds_in (str): The input GeoJSON file path
-        out_conf (dict): The output properties. Must include the keys:
-                         out_file: The output file
-                         size: The output image file (x,y)
-                         geotransform: The output file geotransform
-    """
-    if not all(k in out_conf for k in ("out_file", "size", "geotransform")):
-        raise ValueError("The out_conf parameter doesn't have all" +
-                         "the needed elements")
-
-    layer = ds_in.GetLayer()
-    num_layers = layer.GetFeatureCount()
-    proj = layer.GetSpatialRef()
-
-    driver = gdal.GetDriverByName('GTIFF')
-    ds_out = driver.Create(out_conf['out_file'], out_conf['size'][0],
-                           out_conf['size'][1],
-                           num_layers, gdal.GDT_Float32)
-    ds_out.SetGeoTransform(out_conf['geotransform'])
-    ds_out.SetProjection(proj.ExportToWkt())
-
-    for i in range(num_layers):
-        layer.SetAttributeFilter('cluster={}'.format(i))
-        gdal.RasterizeLayer(ds_out, [i + 1], layer, burn_values=[1])
-
-    data = ds_out.ReadAsArray().astype(numpy.float32)
-    for i in range(num_layers):
-        blurred = gaussian_filter(data[i], sigma)
-        ds_out.GetRasterBand(i + 1).WriteArray(blurred)
-
-    ds_out = None
diff --git a/pymica/utils/create_clusters.py b/pymica/utils/create_clusters.py
deleted file mode 100644
index 8058e78..0000000
--- a/pymica/utils/create_clusters.py
+++ /dev/null
@@ -1,90 +0,0 @@
-'''Tool for creating and editing clusters with the K-Means method
-'''
-import http.server
-import os
-import socketserver
-import webbrowser
-from json import dumps, load
-from math import floor
-from shutil import copytree
-from tempfile import mkdtemp
-
-import numpy as np
-from pyproj import Proj
-from sklearn.cluster import KMeans
-
-
-def create_clusters(locations, n_clusters):
-    """Main function. Creates the clusters and calls the web browser
-       to edit the result
-
-    Args:
-        locations (list, str): The locations list or the file path with the
-                               json data
-        n_clusters (int): The number of clusters to create
-    """
-    if isinstance(locations, str):
-        f_p = open(locations)
-        locations = load(f_p)
-
-    positions = np.zeros([len(locations), 2])
-    utm = calculate_utm_def([locations[0]['lon'], locations[0]['lat']])
-    positions_list = []
-    for i, loc_value in enumerate(locations):
-        positions[i][0] = loc_value['lon']
-        positions[i][1] = loc_value['lat']
-        value = utm(loc_value['lon'], loc_value['lat'])
-        positions_list.append(value)
-    positions_list = np.array(positions_list)
-    kmeans = KMeans(n_clusters=n_clusters, random_state=0).fit(positions_list)
-
-    out_geojson = {"type": "FeatureCollection",
-                   "features": []}
-    for i, labels_value in enumerate(kmeans.labels_):
-        locations[i]['cluster'] = labels_value
-        out_geojson['features'].append({
-            "type": "Feature",
-            "geometry": {
-                "type": "Point",
-                "coordinates": [locations[i]['lon'], locations[i]['lat']]
-            },
-            "properties": {
-                "cluster": int(labels_value),
-                "id": locations[i]['id'],
-                "alt": locations[i]['alt']
-            }
-        })
-
-    web_dir = mkdtemp()
-    copytree("./cluster/web", web_dir+"/web")
-    f_p = open(web_dir+"/web/points.json", "w")
-    f_p.write(dumps(out_geojson))
-    f_p.close()
-
-    os.chdir(web_dir + "/web")
-    socketserver.socket.setdefaulttimeout(5)
-    Handler = http.server.SimpleHTTPRequestHandler
-    httpd = socketserver.TCPServer(("", 8000), Handler)
-    print("serving at port", 8000)
-    webbrowser.get('firefox').open_new_tab("localhost:8000")
-
-    httpd.serve_forever()
-
-
-def calculate_utm_def(point):
-    """Returns the proj4 definition, calculating the utm zone from
-    a point definition
-
-    Args:
-        point (list): A two element list with the lon and lat values
-                      from a point
-
-    Returns:
-        object: The pyproj.Proj object to use to convert
-    """
-
-    zone = (floor((point[0] + 180)/6) % 60) + 1
-    south = " +south " if point[1] < 0 else " "
-    desc = ("+proj=utm +zone={}{}+ellps=WGS84 +datum=WGS84 " +
-            "+units=m +no_defs").format(zone, south)
-    return Proj(desc)
diff --git a/pymica/utils/geotools.py b/pymica/utils/geotools.py
index fd3ac3a..f8bc7c3 100644
--- a/pymica/utils/geotools.py
+++ b/pymica/utils/geotools.py
@@ -1,21 +1,44 @@
+"""Module with geographic tools.
+"""
 import math
+
 import pyproj
 
 
-def calculate_utm_def(point):
-    """Returns the proj4 definition, calculating the utm zone from
-    a point definition
+def get_utm_epsg_from_lonlat(longitude: float, latitude: float) -> str:
+    """Get UTM EPSG from longitude and latitude coordinates.
 
     Args:
-        point (list): A two element list with the lon and lat values
-                      from a point
+        longitude (float): Longitude of a location.
+        latitude (float): Latitude of a location.
 
     Returns:
-        object: The pyproj.Proj object to use to convert
+        str: UTM EPSG coordinates reference system.
     """
+    zone = (math.floor((longitude + 180) / 6) % 60) + 1
+    south = " +south " if latitude < 0 else " "
+    desc = (
+        "+proj=utm +zone={}{}+ellps=WGS84 +datum=WGS84 +units=m +no_defs"
+    ).format(zone, south)
+
+    utm_epsg = pyproj.CRS.from_proj4(desc).to_epsg()
+
+    return utm_epsg
+
+
+def reproject_point(point: tuple, in_proj: int | str, out_proj: int | str) -> tuple:
+    """Reproject a location from `in_proj` to `out_proj`.
+
+    Args:
+        point (tuple): Coordinates of a location.
+        in_proj (int | str): Input projection as EPSG code or proj4.
+        out_proj (int | str): Output projection as EPSG code or proj4.
+
+    Returns:
+        tuple: Reprojected point coordinates.
+    """
+    reproject = pyproj.Transformer.from_crs(in_proj, out_proj, always_xy=True)
+
+    point_x, point_y = reproject.transform(point[0], point[1])
 
-    zone = (math.floor((point[0] + 180)/6) % 60) + 1
-    south = " +south " if point[1] < 0 else " "
-    desc = ("+proj=utm +zone={}{}+ellps=WGS84 +datum=WGS84 " +
-            "+units=m +no_defs").format(zone, south)
-    return pyproj.Proj(desc)
\ No newline at end of file
+    return (point_x, point_y)


Test patch:
diff --git a/pymica_tests/clusters_test.py b/pymica_tests/clusters_test.py
new file mode 100644
index 0000000..db653db
--- /dev/null
+++ b/pymica_tests/clusters_test.py
@@ -0,0 +1,165 @@
+"""Tests creation of clusters files (shapefiles and rasters)
+"""
+import unittest
+from os import remove
+from os.path import exists
+from tempfile import gettempdir
+
+from osgeo import gdal, ogr
+import json
+
+from pymica.utils.clusters import (
+    create_reprojected_geometries,
+    rasterize_clusters,
+    create_clusters,
+)
+
+
+class TestClusters(unittest.TestCase):
+    """Test creation and rasterization of clusters"""
+
+    with open("pymica_tests/data/sample_station_metadata.json", "r") as f_p:
+        sample_data = json.load(f_p)
+
+    def test_create_reprojected_geometries(self):
+        """Test creation of reprojected geometries"""
+        d_s = create_reprojected_geometries("pymica_tests/data/clusters.json", 25831)
+
+        self.assertIsInstance(d_s, ogr.DataSource)
+        layer = d_s.GetLayer()
+
+        self.assertEqual(layer.GetFeatureCount(), 3)
+        self.assertTrue(layer.GetSpatialRef() is not None)
+
+        with self.assertRaises(IOError) as err:
+            create_reprojected_geometries("/bad_file.json", 25831)
+        self.assertEqual("File /bad_file.json doesn't exist", str(err.exception))
+
+        with self.assertRaises(ValueError) as err:
+            create_reprojected_geometries("pymica_tests/data/clusters.json", -999)
+        self.assertEqual("Wrong EPSG code: -999", str(err.exception))
+
+    def test_rasterize_clusters(self):
+        """Test rasterization of clusters"""
+        out_file = gettempdir() + "/test.tiff"
+        geotransform = (260000, 270, 0, 4750000, 0, -270)
+        if exists(out_file):
+            remove(out_file)
+
+        layer = create_reprojected_geometries("pymica_tests/data/clusters.json", 25831)
+        rasterize_clusters(
+            layer,
+            {"size": [1000, 1000], "geotransform": geotransform, "out_file": out_file},
+        )
+        d_s = gdal.Open(out_file)
+        self.assertIsNotNone(d_s)
+
+        self.assertEqual(d_s.GetGeoTransform(), geotransform)
+        self.assertEqual(d_s.RasterCount, 3)
+
+        # Errors
+        with self.assertRaises(KeyError) as err:
+            rasterize_clusters(
+                layer,
+                {
+                    "size": [1000, 1000],
+                    "geotransform": geotransform,
+                },
+            )
+        self.assertEqual(
+            "`out_file`, `size`, `geotransform` must be in the `raster_config` "
+            "parameter.",
+            err.exception.args[0],
+        )
+
+        with self.assertRaises(KeyError) as err:
+            rasterize_clusters(
+                layer,
+                {"out_file": gettempdir() + "/test.tiff", "geotransform": geotransform},
+            )
+        self.assertEqual(
+            "`out_file`, `size`, `geotransform` must be in the `raster_config` "
+            "parameter.",
+            err.exception.args[0],
+        )
+
+        with self.assertRaises(KeyError) as err:
+            rasterize_clusters(
+                layer, {"out_file": gettempdir() + "/test.tiff", "size": [1000, 1000]}
+            )
+        self.assertEqual(
+            "`out_file`, `size`, `geotransform` must be in the `raster_config` "
+            "parameter.",
+            err.exception.args[0],
+        )
+
+    def test_create_clusters(self):
+        """Test creation of clusters with Kmeans"""
+        create_clusters(self.sample_data, 6, "pymica_tests/data/test_clusters_6.json")
+
+        self.assertTrue(exists("pymica_tests/data/test_clusters_6.json"))
+
+        with open("pymica_tests/data/test_clusters_6.json", "r") as f_p:
+            clusters = json.load(f_p)
+
+        self.assertEqual(
+            clusters["features"][0]["geometry"]["coordinates"], [2.18091, 41.39004]
+        )
+        self.assertEqual(
+            clusters["features"][0]["properties"],
+            {"cluster": 1, "id": "AN", "alt": 7.5},
+        )
+
+        self.assertEqual(
+            clusters["features"][-1]["geometry"]["coordinates"], [1.89716, 42.32211]
+        )
+        self.assertEqual(
+            clusters["features"][-1]["properties"],
+            {"cluster": 2, "id": "ZD", "alt": 2478},
+        )
+
+        if exists("pymica_tests/data/test_clusters_6.json"):
+            remove("pymica_tests/data/test_clusters_6.json")
+
+    def test_create_clusters_from_file(self):
+        """Test creation of clusters with Kmeans from file"""
+        create_clusters(
+            "pymica_tests/data/sample_station_metadata.json",
+            6,
+            "pymica_tests/data/test_clusters_6.json",
+        )
+
+        self.assertTrue(exists("pymica_tests/data/test_clusters_6.json"))
+
+        with open("pymica_tests/data/test_clusters_6.json", "r") as f_p:
+            clusters = json.load(f_p)
+
+        self.assertEqual(
+            clusters["features"][0]["geometry"]["coordinates"], [2.18091, 41.39004]
+        )
+        self.assertEqual(
+            clusters["features"][0]["properties"],
+            {"cluster": 1, "id": "AN", "alt": 7.5},
+        )
+
+        self.assertEqual(
+            clusters["features"][-1]["geometry"]["coordinates"], [1.89716, 42.32211]
+        )
+        self.assertEqual(
+            clusters["features"][-1]["properties"],
+            {"cluster": 2, "id": "ZD", "alt": 2478},
+        )
+
+        if exists("pymica_tests/data/test_clusters_6.json"):
+            remove("pymica_tests/data/test_clusters_6.json")
+
+    def test_create_clusters_wrong_dir(self):
+        """Test error raise from creation of clusters"""
+        with self.assertRaises(FileNotFoundError) as err:
+            create_clusters(
+                "pymica_tests/data/sample_station_metadata.json",
+                6,
+                "not/found/test_clusters_6.json",
+            )
+
+        self.assertEqual(err.exception.args[0], "not/found directory does not exist.")
diff --git a/pymica_tests/create_cluster_files_test.py b/pymica_tests/create_cluster_files_test.py
deleted file mode 100644
index 07fd3b8..0000000
--- a/pymica_tests/create_cluster_files_test.py
+++ /dev/null
@@ -1,79 +0,0 @@
-'''Tests fot create_cluster_files.py
-'''
-import unittest
-from os import remove
-from os.path import exists
-from tempfile import gettempdir
-
-from osgeo import gdal, ogr
-
-from pymica.utils.create_cluster_files import (create_reprojected_geoms,
-                                          rasterize_clusters)
-
-
-class TestCreateClustersFiles(unittest.TestCase):
-    def test_create_reprojected_geoms(self):
-        d_s = create_reprojected_geoms("pymica_tests/data/clusters.json", 25831)
-
-        self.assertIsInstance(d_s, ogr.DataSource)
-        layer = d_s.GetLayer()
-
-        self.assertEqual(layer.GetFeatureCount(), 3)
-        self.assertTrue(layer.GetSpatialRef() is not None)
-
-        with self.assertRaises(IOError) as c_m:
-            create_reprojected_geoms("/bad_file.json", 25831)
-        self.assertEqual(
-            "File /bad_file.json doesn't exist",
-            str(c_m.exception))
-
-        with self.assertRaises(ValueError) as c_m:
-            create_reprojected_geoms("pymica_tests/data/clusters.json", -999)
-        self.assertEqual(
-            "Wrong EPSG code: -999",
-            str(c_m.exception))
-
-    def test_rasterize_clusters(self):
-        out_file = gettempdir() + "/test.tiff"
-        geotransform = (260000, 270, 0, 4750000, 0, -270)
-        if exists(out_file):
-            remove(out_file)
-
-        layer = create_reprojected_geoms("pymica_tests/data/clusters.json", 25831)
-        rasterize_clusters(layer,
-                           {'size': [1000, 1000],
-                            'geotransform': geotransform,
-                            'out_file': out_file})
-        d_s = gdal.Open(out_file)
-        self.assertIsNotNone(d_s)
-
-        self.assertEqual(d_s.GetGeoTransform(), geotransform)
-        self.assertEqual(d_s.RasterCount, 3)
-
-        # Errors
-        with self.assertRaises(ValueError) as c_m:
-            rasterize_clusters(layer,
-                               {'size': [1000, 1000],
-                                'geotransform': geotransform,
-                                })
-        self.assertEqual(
-            "The out_conf parameter doesn't have allthe needed elements",
-            str(c_m.exception))
-
-        with self.assertRaises(ValueError) as c_m:
-            rasterize_clusters(layer,
-                               {'out_file': gettempdir() + '/test.tiff',
-                                'geotransform': geotransform
-                                })
-        self.assertEqual(
-            "The out_conf parameter doesn't have allthe needed elements",
-            str(c_m.exception))
-
-        with self.assertRaises(ValueError) as c_m:
-            rasterize_clusters(layer,
-                               {'out_file': gettempdir() + '/test.tiff',
-                                'size': [1000, 1000]
-                                })
-        self.assertEqual(
-            "The out_conf parameter doesn't have allthe needed elements",
-            str(c_m.exception))
diff --git a/pymica_tests/create_clusters_test.py b/pymica_tests/create_clusters_test.py
deleted file mode 100644
index 66fdc82..0000000
--- a/pymica_tests/create_clusters_test.py
+++ /dev/null
@@ -1,46 +0,0 @@
-'''Tests for cluster.create_clusters.py
-'''
-import unittest
-
-from pymica.utils.create_clusters import calculate_utm_def
-from pyproj import Proj
-
-
-class TestCreateClusters(unittest.TestCase):
-
-    # Not possible to open a browser in TravisCI
-    '''
-    def test_create_clusters(self):
-        fp = open("./test/data/sample_station_metadata.json")
-        data = json.load(fp)
-        create_clusters(data, 6)
-        fp.close()
-    '''
-
-    def test_calculate_utm_def(self):
-        point = [0, 0]
-        result = calculate_utm_def(point)
-        expected = Proj("+proj=utm +zone=31 +ellps=WGS84 " +
-                        "+datum=WGS84 +units=m +no_defs")
-
-        self.assertIsInstance(result, Proj)
-        self.assertEqual(set(result(point[0], point[1])),
-                         set(expected(point[0], point[1])))
-
-        point = [0, -10]
-        result = calculate_utm_def(point)
-        expected = Proj("+proj=utm +zone=31 +south +ellps=WGS84 " +
-                        "+datum=WGS84 +units=m +no_defs")
-
-        self.assertIsInstance(result, Proj)
-        self.assertEqual(set(result(point[0], point[1])),
-                         set(expected(point[0], point[1])))
-
-        point = [-178, 24]
-        result = calculate_utm_def(point)
-        expected = Proj("+proj=utm +zone=1 +ellps=WGS84 " +
-                        "+datum=WGS84 +units=m +no_defs")
-
-        self.assertIsInstance(result, Proj)
-        self.assertEqual(set(result(point[0], point[1])),
-                         set(expected(point[0], point[1])))
diff --git a/pymica_tests/geotools_test.py b/pymica_tests/geotools_test.py
new file mode 100644
index 0000000..7cedc19
--- /dev/null
+++ b/pymica_tests/geotools_test.py
@@ -0,0 +1,30 @@
+"""Tests for creation of clusters.
+"""
+import unittest
+
+import pyproj
+
+from pymica.utils.geotools import get_utm_epsg_from_lonlat
+
+
+class TestGeotools(unittest.TestCase):
+    """Test geotools module"""
+    def test_get_utm_epsg_from_lonlat(self):
+        """Test get definition of UTM from a point"""
+        utm_epsg = get_utm_epsg_from_lonlat(0, 0)
+        expected = pyproj.CRS.from_proj4(
+            "+proj=utm +zone=31 +ellps=WGS84 +datum=WGS84 +units=m +no_defs"
+        ).to_epsg()
+        self.assertEqual(expected, utm_epsg)
+
+        utm_epsg = get_utm_epsg_from_lonlat(0, -10)
+        expected = pyproj.CRS.from_proj4(
+            "+proj=utm +zone=31 +south +ellps=WGS84 +datum=WGS84 +units=m +no_defs"
+        ).to_epsg()
+        self.assertEqual(expected, utm_epsg)
+
+        utm_epsg = get_utm_epsg_from_lonlat(-178, 24)
+        expected = pyproj.CRS.from_proj4(
+            "+proj=utm +zone=1 +ellps=WGS84 +datum=WGS84 +units=m +no_defs"
+        ).to_epsg()
+        self.assertEqual(expected, utm_epsg)