Support .gpkg files in shp options

contribs/application/src/main/java/org/matsim/application/options/ShpOptions.java

@@ -1,11 +1,12 @@
 package org.matsim.application.options;
 
+import org.apache.commons.lang3.ArrayUtils;
 import org.apache.logging.log4j.LogManager;
 import org.apache.logging.log4j.Logger;
-import org.geotools.data.FeatureReader;
-import org.geotools.data.FileDataStoreFactorySpi;
+import org.geotools.data.*;
 import org.geotools.data.shapefile.ShapefileDataStore;
 import org.geotools.data.shapefile.ShapefileDataStoreFactory;
+import org.geotools.geopkg.GeoPkgDataStoreFactory;
 import org.geotools.referencing.CRS;
 import org.locationtech.jts.geom.*;
 import org.locationtech.jts.index.strtree.AbstractNode;
@@ -34,10 +35,7 @@
 import java.net.URL;
 import java.nio.charset.Charset;
 import java.nio.file.Path;
-import java.util.ArrayList;
-import java.util.Collection;
-import java.util.List;
-import java.util.Set;
+import java.util.*;
 import java.util.function.Predicate;
 import java.util.zip.ZipEntry;
 import java.util.zip.ZipInputStream;
@@ -62,6 +60,9 @@ public final class ShpOptions {
 	@CommandLine.Option(names = "--shp-crs", description = "Overwrite coordinate system of the shape file")
 	private String shpCrs;
 
+	@CommandLine.Option(names = "--shp-layer", description = "Layer to use (for .gpgk files). Defaults to first layer.", required = false)
+	private String shpLayer;
+
 	@CommandLine.Option(names = "--shp-charset", description = "Charset used to read the shape file", defaultValue = "ISO-8859-1")
 	private Charset shpCharset;
 
@@ -86,19 +87,41 @@ public ShpOptions(String shp, @Nullable String shpCrs, @Nullable Charset shpChar
 		this.shpCharset = shpCharset;
 	}
 
+	private ShpOptions(String shp, String shpCrs, String shpLayer, Charset shpCharset) {
+		this.shp = shp;
+		this.shpCrs = shpCrs;
+		this.shpLayer = shpLayer;
+		this.shpCharset = shpCharset;
+	}
+
+	/**
+	 * Create shp options with a specific layer. (Usually for gpkg files).
+	 */
+	public static ShpOptions ofLayer(String shp, @Nullable String shpLayer) {
+		return new ShpOptions(shp, null, shpLayer, null);
+	}
+
 	/**
 	 * Opens datastore to a shape-file.
 	 */
-	public static ShapefileDataStore openDataStore(String shp) throws IOException {
+	public static DataStore openDataStore(String shp) throws IOException {
 
 		FileDataStoreFactorySpi factory = new ShapefileDataStoreFactory();
 
 		URL url = IOUtils.resolveFileOrResource(shp);
 
-		ShapefileDataStore ds;
+		DataStore ds;
 		if (shp.endsWith(".shp"))
-			ds = (ShapefileDataStore) factory.createDataStore(url);
-		else if (shp.endsWith(".zip")) {
+			ds = factory.createDataStore(url);
+		else if (shp.endsWith(".gpkg")) {
+			ds = DataStoreFinder.getDataStore(Map.of(
+				GeoPkgDataStoreFactory.DBTYPE.key, "geopkg",
+				GeoPkgDataStoreFactory.DATABASE.key, shp,
+				GeoPkgDataStoreFactory.READ_ONLY.key, true
+			));
+		} else if (shp.endsWith(".dbf") || shp.endsWith(".shx")) {
+			throw new IllegalArgumentException("Shape file must be .shp, but was: " + shp);
+		} else if (shp.endsWith(".zip")) {
 
 			// Zip files will only work with local files
 			URI uri;
@@ -119,12 +142,12 @@ else if (shp.endsWith(".zip")) {
 				throw new IllegalArgumentException("Could not create URI for zip file: " + url, e);
 			}
 
-			ds = (ShapefileDataStore) factory.createDataStore(uri.toURL());
+			ds = factory.createDataStore(uri.toURL());
 		} else {
 			throw new IllegalArgumentException("Shape file must either be .zip or .shp, but was: " + shp);
 		}
 
-		return ds;
+		return Objects.requireNonNull(ds, "Could not create data store.");
 	}
 
 	/**
@@ -158,11 +181,15 @@ public List<SimpleFeature> readFeatures() {
 			throw new IllegalStateException("Shape file path not specified");
 
 		try {
-			ShapefileDataStore ds = openDataStore(shp);
-			if (shpCharset != null)
-				ds.setCharset(shpCharset);
+			DataStore ds = openDataStore(shp);
+			if (shpCharset != null && ds instanceof ShapefileDataStore shpDs)
+				shpDs.setCharset(shpCharset);
 
-			return GeoFileReader.getSimpleFeatures(ds);
+			if (ds instanceof FileDataStore fds) {
+				return GeoFileReader.getSimpleFeatures(fds);
+			}
+
+			return GeoFileReader.getSimpleFeatures(ds, getLayer(ds));
 		} catch (IOException e) {
 			throw new UncheckedIOException(e);
 		}
@@ -179,7 +206,7 @@ public Geometry getGeometry() {
 		}
 
 		if (features.size() == 1) {
-            return (Geometry) features.iterator().next().getDefaultGeometry();
+			return (Geometry) features.iterator().next().getDefaultGeometry();
 		}
 
 		GeometryFactory factory = ((Geometry) features.iterator().next().getDefaultGeometry()).getFactory();
@@ -212,7 +239,7 @@ public Index createIndex(String queryCRS, String attr, Predicate<SimpleFeature>
 			throw new IllegalArgumentException("Query crs must not be null!");
 
 		try {
-			ShapefileDataStore ds = openDataStoreAndSetCRS();
+			DataStore ds = openDataStoreAndSetCRS();
 			CoordinateTransformation ct;
 			if (queryCRS.equals(SAME_CRS))
 				ct = new IdentityTransformation();
@@ -254,8 +281,8 @@ private String detectCRS() {
 
 		if (shpCrs == null) {
 			try {
-				ShapefileDataStore ds = openDataStore(shp);
-				CoordinateReferenceSystem crs = ds.getSchema().getCoordinateReferenceSystem();
+				DataStore ds = openDataStore(shp);
+				CoordinateReferenceSystem crs = ds.getSchema(getLayer(ds)).getCoordinateReferenceSystem();
 				ds.dispose();
 				shpCrs = "EPSG:" + CRS.lookupEpsgCode(crs, true);
 				log.info("Using detected crs for {}: {}", shp, shpCrs);
@@ -271,12 +298,12 @@ private String detectCRS() {
 	/**
 	 * Open the shape file for processing and set the crs if not already specified.
 	 */
-	private ShapefileDataStore openDataStoreAndSetCRS() throws IOException {
-		ShapefileDataStore ds = openDataStore(shp);
+	private DataStore openDataStoreAndSetCRS() throws IOException {
+		DataStore ds = openDataStore(shp);
 
 		if (shpCrs == null) {
 			try {
-				CoordinateReferenceSystem crs = ds.getSchema().getCoordinateReferenceSystem();
+				CoordinateReferenceSystem crs = ds.getSchema(getLayer(ds)).getCoordinateReferenceSystem();
 				shpCrs = "EPSG:" + CRS.lookupEpsgCode(crs, true);
 				log.info("Using detected crs for {}: {}", shp, shpCrs);
 			} catch (FactoryException | NullPointerException e) {
@@ -287,6 +314,21 @@ private ShapefileDataStore openDataStoreAndSetCRS() throws IOException {
 		return ds;
 	}
 
+	/**
+	 * Get the selected layer or throw exception if not found.
+	 */
+	private String getLayer(DataStore ds) throws IOException {
+		String[] typeNames = ds.getTypeNames();
+		if (shpLayer != null) {
+			if (!ArrayUtils.contains(typeNames, shpLayer))
+				throw new IllegalArgumentException("Layer " + shpLayer + " not found in shape file.");
+
+			return shpLayer;
+		}
+
+		return typeNames[0];
+	}
+
 	/**
 	 * Create an inverse coordinate transformation from the shape file crs. Tries to autodetect the crs of the shape file.
 	 */
@@ -309,11 +351,18 @@ public final class Index {
 		 * @param ct   coordinate transform from query to target crs
 		 * @param attr attribute for the result of {@link #query(Coord)}
 		 */
-		Index(CoordinateTransformation ct, ShapefileDataStore ds, String attr, @Nullable Predicate<SimpleFeature> filter) throws IOException {
-			if (shpCharset != null)
-				ds.setCharset(shpCharset);
+		Index(CoordinateTransformation ct, DataStore ds, String attr, @Nullable Predicate<SimpleFeature> filter) throws IOException {
+			if (shpCharset != null && ds instanceof ShapefileDataStore shpDs)
+				shpDs.setCharset(shpCharset);
+
+			FeatureReader<SimpleFeatureType, SimpleFeature> it;
+			DefaultTransaction transaction = new DefaultTransaction();
+			if (ds instanceof FileDataStore fds) {
+				it = fds.getFeatureReader();
+			} else {
+				it = ds.getFeatureReader(new Query(getLayer(ds)), transaction);
+			}
 
-			FeatureReader<SimpleFeatureType, SimpleFeature> it = ds.getFeatureReader();
 			while (it.hasNext()) {
 				SimpleFeature ft = it.next();
 
@@ -332,6 +381,7 @@ public final class Index {
 
 			index.build();
 
+			transaction.close();
 			it.close();
 			ds.dispose();
 

contribs/application/src/main/java/org/matsim/application/prepare/Predictor.java

@@ -0,0 +1,32 @@
+package org.matsim.application.prepare;
+
+import it.unimi.dsi.fastutil.objects.Object2DoubleMap;
+import it.unimi.dsi.fastutil.objects.Object2ObjectMap;
+
+/**
+ * Predictor interface for arbitrary numeric values.
+ */
+public interface Predictor {
+
+
+	/**
+	 * Predict value from given features.
+	 */
+	double predict(Object2DoubleMap<String> features, Object2ObjectMap<String, String> categories);
+
+	/**
+	 * Predict values with adjusted model params.
+	 */
+	default double predict(Object2DoubleMap<String> features, Object2ObjectMap<String, String> categories, double[] params) {
+		throw new UnsupportedOperationException("Not implemented");
+	}
+
+
+	/**
+	 * Return data that is used for internal prediction function (normalization already applied).
+	 */
+	default double[] getData(Object2DoubleMap<String> features, Object2ObjectMap<String, String> categories) {
+		throw new UnsupportedOperationException("Not implemented");
+	}
+
+}

contribs/application/src/test/java/org/matsim/application/options/ShpOptionsTest.java

@@ -29,7 +29,7 @@ void readZip() {
 		Path input = Path.of(utils.getClassInputDirectory()
 				.replace("ShpOptionsTest", "CreateLandUseShpTest")
 				.replace("options", "prepare"))
-				.resolve("andorra-latest-free.shp.zip");
+			.resolve("andorra-latest-free.shp.zip");
 
 		Assumptions.assumeTrue(Files.exists(input));
 
@@ -38,7 +38,7 @@ void readZip() {
 		List<SimpleFeature> ft = shp.readFeatures();
 
 		assertThat(ft)
-				.isNotEmpty();
+			.isNotEmpty();
 
 	}
 
@@ -73,7 +73,7 @@ void all() {
 		Path input = Path.of(utils.getClassInputDirectory()
 				.replace("ShpOptionsTest", "CreateLandUseShpTest")
 				.replace("options", "prepare"))
-				.resolve("andorra-latest-free.shp.zip");
+			.resolve("andorra-latest-free.shp.zip");
 
 		Assumptions.assumeTrue(Files.exists(input));
 
@@ -84,8 +84,8 @@ void all() {
 		List<SimpleFeature> ft = index.getAllFeatures();
 
 		assertThat(ft)
-				.hasSize(4906)
-				.hasSize(Set.copyOf(ft).size());
+			.hasSize(4906)
+			.hasSize(Set.copyOf(ft).size());
 
 		assertThat(shp.readFeatures())
 			.hasSameElementsAs(ft);
@@ -96,17 +96,36 @@ void all() {
 	void testGetGeometry() {
 
 		Path input = Path.of(utils.getClassInputDirectory()
-						.replace("ShpOptionsTest", "CreateLandUseShpTest")
-						.replace("options", "prepare"))
-				.resolve("andorra-latest-free.shp.zip");
+				.replace("ShpOptionsTest", "CreateLandUseShpTest")
+				.replace("options", "prepare"))
+			.resolve("andorra-latest-free.shp.zip");
 
 		Assumptions.assumeTrue(Files.exists(input));
 
 		ShpOptions shp = new ShpOptions(input, null, null);
-		Geometry geometry = shp.getGeometry() ;
+		Geometry geometry = shp.getGeometry();
 
 		assertThat(geometry.getArea())
 			.isCloseTo(1.9847543618489646E-4, Offset.offset(1e-8));
 
 	}
+
+	@Test
+	void gpkg() {
+
+		Path path = Path.of(utils.getPackageInputDirectory(), "example.gpkg");
+
+		ShpOptions shp = ShpOptions.ofLayer(path.toString(), null);
+
+		List<SimpleFeature> features = shp.readFeatures();
+
+		assertThat(features)
+			.hasSize(3);
+
+		ShpOptions.Index index = shp.createIndex("_");
+
+		assertThat(index.size())
+			.isEqualTo(3);
+
+	}
 }

matsim/src/main/java/org/matsim/core/utils/gis/GeoFileReader.java

@@ -133,6 +133,19 @@ public static List<SimpleFeature> getSimpleFeatures(DataStore dataStore, Name la
 		return featureSet;
 	}
 
+	/**
+	 * Read all simple features from a data store. This method makes sure the store is closed afterwards.
+	 * @return list of contained features.
+	 * @see #getSimpleFeatures(DataStore, Name)
+	 */
+	public static List<SimpleFeature> getSimpleFeatures(DataStore dataStore, String layerName) throws IOException {
+		SimpleFeatureSource featureSource = dataStore.getFeatureSource(layerName);
+		Gbl.assertNotNull(featureSource);
+		List<SimpleFeature> featureSet = getSimpleFeatures(featureSource);
+		dataStore.dispose();
+		return featureSet;
+	}
+
 	private static List<SimpleFeature> getSimpleFeatures(SimpleFeatureSource featureSource) throws IOException {
 		SimpleFeatureIterator it = featureSource.getFeatures().features();
 		List<SimpleFeature> featureSet = new ArrayList<>();