Geotools 拆分光栅并再次合并它不起作用
Geotools splitting raster and merging it again doesn't work
第一次在这里发帖，如有错误请见谅。
我们正在尝试让用户选择一张光栅图(.tiff)，以便上传到hbase进行处理等操作。由于单张图片可能很大，我们尝试将它们拆分，以便加快上传过程，然后再将它们合并在一起，得到原始的单张图片。
据我所知明白了，错误似乎出在split方法上。
目前我们只是尝试拆分镜像并在本地重建它。
这是拆分方法的代码
public ArrayList<GridCoverage2D> createTiles(GridCoverage2D grid) throws FactoryException, TransformException {
Geometry geom = buildGeometryFromGrid(grid);
org.locationtech.jts.geom.Envelope env = geom.getEnvelopeInternal();
Coverage coverage = GeoHash.coverBoundingBoxMaxHashes(env.getMaxY(), env.getMinX(), env.getMinY(), env.getMaxX(), 1);
ArrayList<String> geohashes = new ArrayList<>();
String hash = null;
if (coverage == null) hash = "";
else {
for (String singleHash : coverage.getHashes()) {
hash = singleHash;
break;
}
}
if (hash.length() == Constants.defaultPrecision) {
coverage = GeoHash.coverBoundingBox(env.getMaxY(), env.getMinX(), env.getMinY(), env.getMaxX(), Constants.defaultPrecision + 1);
geohashes = new ArrayList<>(coverage.getHashes());
}
else if (hash.length() < Constants.defaultPrecision) {
coverage = GeoHash.coverBoundingBox(env.getMaxY(), env.getMinX(), env.getMinY(), env.getMaxX(), Constants.defaultPrecision);
geohashes = new ArrayList<>(coverage.getHashes());
}
else {
if (hash.length() > Constants.spaceLength) geohashes.add(hash.substring(0, hash.length() - 1));
else if (hash.length() == Constants.spaceLength) geohashes.add(hash);
else {
coverage = GeoHash.coverBoundingBox(env.getMaxY(), env.getMinX(), env.getMinY(), env.getMaxX(), hash.length() + 1);
geohashes = new ArrayList<>(coverage.getHashes());
}
}
System.out.println(geohashes.size());
List<ReferencedEnvelope> envelopes = calcuateCuts(geohashes);
ArrayList<GridCoverage2D> tiles = new ArrayList<>();
for (ReferencedEnvelope envelope : envelopes) {
tiles.addAll(getTiles(envelope, grid));
}
System.out.println(tiles.size());
return tiles;
}
它有点管用，因为我们得到了一些 tiff，但我们也得到了其中一些无法使用 qgis 等工具可视化的问题。 merge 方法给了我们一个 tiff 图像，但是它又没有在 qgis.
中可视化
如果您需要更多信息，请告诉我。
谢谢！
编辑
这是原始且正确的 tif
的 gdalinfo
Driver: GTiff/GeoTIFF
Files: C:\Users\Consulente\Desktop\appoggio_rilascio\originale.tif
C:\Users\Consulente\Desktop\appoggio_rilascio\originale.tif.aux.xml
Size is 266, 269
Coordinate System is:
PROJCRS["WGS 84 / UTM zone 33N",
BASEGEOGCRS["WGS 84",
DATUM["World Geodetic System 1984",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4326]],
CONVERSION["UTM zone 33N",
METHOD["Transverse Mercator",
ID["EPSG",9807]],
PARAMETER["Latitude of natural origin",0,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8801]],
PARAMETER["Longitude of natural origin",15,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Scale factor at natural origin",0.9996,
SCALEUNIT["unity",1],
ID["EPSG",8805]],
PARAMETER["False easting",500000,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",0,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["(E)",east,
ORDER[1],
LENGTHUNIT["metre",1]],
AXIS["(N)",north,
ORDER[2],
LENGTHUNIT["metre",1]],
USAGE[
SCOPE["unknown"],
AREA["World - N hemisphere - 12┬░E to 18┬░E - by country"],
BBOX[0,12,84,18]],
ID["EPSG",32633]]
Data axis to CRS axis mapping: 1,2
Origin = (293436.000000000000000,4647354.000000000000000)
Pixel Size = (3.000000000000000,-3.000000000000000)
Metadata:
AREA_OR_POINT=Area
Image Structure Metadata:
COMPRESSION=LZW
INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left ( 293436.000, 4647354.000) ( 12d30'28.04"E, 41d57' 4.05"N)
Lower Left ( 293436.000, 4646547.000) ( 12d30'29.05"E, 41d56'37.91"N)
Upper Right ( 294234.000, 4647354.000) ( 12d31' 2.66"E, 41d57' 4.80"N)
Lower Right ( 294234.000, 4646547.000) ( 12d31' 3.68"E, 41d56'38.66"N)
Center ( 293835.000, 4646950.500) ( 12d30'45.86"E, 41d56'51.36"N)
Band 1 Block=266x3 Type=UInt16, ColorInterp=Red
Min=4230.000 Max=18752.000
Minimum=4230.000, Maximum=18752.000, Mean=6493.946, StdDev=1353.641
NoData Value=0
Metadata:
STATISTICS_MAXIMUM=18752
STATISTICS_MEAN=6493.9456647057
STATISTICS_MINIMUM=4230
STATISTICS_STDDEV=1353.6413052408
STATISTICS_VALID_PERCENT=94.55
Band 2 Block=266x3 Type=UInt16, ColorInterp=Green
Min=3979.000 Max=17611.000
Minimum=3979.000, Maximum=17611.000, Mean=5971.019, StdDev=1344.279
NoData Value=0
Metadata:
STATISTICS_MAXIMUM=17611
STATISTICS_MEAN=5971.0194223549
STATISTICS_MINIMUM=3979
STATISTICS_STDDEV=1344.2792863221
STATISTICS_VALID_PERCENT=94.55
Band 3 Block=266x3 Type=UInt16, ColorInterp=Blue
Min=2619.000 Max=16659.000
Minimum=2619.000, Maximum=16659.000, Mean=4982.758, StdDev=1525.882
NoData Value=0
Metadata:
STATISTICS_MAXIMUM=16659
STATISTICS_MEAN=4982.7577379017
STATISTICS_MINIMUM=2619
STATISTICS_STDDEV=1525.8817457971
STATISTICS_VALID_PERCENT=94.55
Band 4 Block=266x3 Type=UInt16, ColorInterp=Undefined
NoData Value=0
这是来自不正确合并 tiff 的 gdalinfo（应该与原始 tiff 相同）
Driver: GTiff/GeoTIFF
Files: C:\Users\Consulente\Desktop\appoggio_rilascio\merged_tif\merged.tiff
Size is 266, 269
Coordinate System is:
GEOGCRS["WGS 84",
DATUM["World Geodetic System 1984",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
CS[ellipsoidal,2],
AXIS["geodetic latitude (Lat)",north,
ORDER[1],
ANGLEUNIT["degree",0.0174532925199433]],
AXIS["geodetic longitude (Lon)",east,
ORDER[2],
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4326]]
Data axis to CRS axis mapping: 2,1
GeoTransform =
41.95133415197179, 0, -2.77698308382438e-05
12.50778763370388, 3.722213354286106e-05, 0
Metadata:
AREA_OR_POINT=Area
TIFFTAG_RESOLUTIONUNIT=1 (unitless)
TIFFTAG_XRESOLUTION=1
TIFFTAG_YRESOLUTION=1
Image Structure Metadata:
INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left ( 41.9513342, 12.5077876) ( 41d57' 4.80"E, 12d30'28.04"N)
Lower Left ( 41.9438641, 12.5077876) ( 41d56'37.91"E, 12d30'28.04"N)
Upper Right ( 41.9513342, 12.5176887) ( 41d57' 4.80"E, 12d31' 3.68"N)
Lower Right ( 41.9438641, 12.5176887) ( 41d56'37.91"E, 12d31' 3.68"N)
Center ( 41.9475991, 12.5127382) ( 41d56'51.36"E, 12d30'45.86"N)
Band 1 Block=266x8 Type=UInt16, ColorInterp=Red
NoData Value=0
Band 2 Block=266x8 Type=UInt16, ColorInterp=Green
NoData Value=0
Band 3 Block=266x8 Type=UInt16, ColorInterp=Blue
NoData Value=0
Band 4 Block=266x8 Type=UInt16, ColorInterp=Alpha
NoData Value=0
通过查看 GdalInfo 输出，您已将 CRS 更改为 EPSG:4326 并翻转了 lat/lon 值 - 所以在您写出您没有的图块的方式中发生了一些事情发布代码。
如果我尝试这样做，我会查看 Image Tiling tutorial，关键的内部循环是这样的：
// iterate over our tile counts
for (int i = 0; i < htc; i++) {
for (int j = 0; j < vtc; j++) {
System.out.println("Processing tile at indices i: " + i + " and j: " + j);
// create the envelope of the tile
Envelope envelope = getTileEnvelope(coverageMinX, coverageMinY, geographicTileWidth, geographicTileHeight,
targetCRS, i, j);
GridCoverage2D finalCoverage = cropCoverage(gridCoverage, envelope);
if (this.getTileScale() != null) {
finalCoverage = scaleCoverage(finalCoverage);
}
// use the AbstractGridFormat's writer to write out the tile
fileExtension = "png";
File tileFile = new File(tileDirectory, i + "_" + j + "." + fileExtension);
final WorldImageWriter wiWriter = new WorldImageWriter(tileFile);
// writing parameters for png
final Format oFormat = wiWriter.getFormat();
((AbstractGridFormat) oFormat).getWriter(tileFile).write(finalCoverage, null);
}
}
和 cropCoverage 使用 CoverageProcessor class:
实现
private GridCoverage2D cropCoverage(GridCoverage2D gridCoverage, Envelope envelope) {
CoverageProcessor processor = CoverageProcessor.getInstance();
// An example of manually creating the operation and parameters we want
final ParameterValueGroup param = processor.getOperation("CoverageCrop").getParameters();
param.parameter("Source").setValue(gridCoverage);
param.parameter("Envelope").setValue(envelope);
return (GridCoverage2D) processor.doOperation(param);
}
第一次在这里发帖，如有错误请见谅。
我们正在尝试让用户选择一张光栅图(.tiff)，以便上传到hbase进行处理等操作。由于单张图片可能很大，我们尝试将它们拆分，以便加快上传过程，然后再将它们合并在一起，得到原始的单张图片。
据我所知明白了，错误似乎出在split方法上。
目前我们只是尝试拆分镜像并在本地重建它。
这是拆分方法的代码
public ArrayList<GridCoverage2D> createTiles(GridCoverage2D grid) throws FactoryException, TransformException {
Geometry geom = buildGeometryFromGrid(grid);
org.locationtech.jts.geom.Envelope env = geom.getEnvelopeInternal();
Coverage coverage = GeoHash.coverBoundingBoxMaxHashes(env.getMaxY(), env.getMinX(), env.getMinY(), env.getMaxX(), 1);
ArrayList<String> geohashes = new ArrayList<>();
String hash = null;
if (coverage == null) hash = "";
else {
for (String singleHash : coverage.getHashes()) {
hash = singleHash;
break;
}
}
if (hash.length() == Constants.defaultPrecision) {
coverage = GeoHash.coverBoundingBox(env.getMaxY(), env.getMinX(), env.getMinY(), env.getMaxX(), Constants.defaultPrecision + 1);
geohashes = new ArrayList<>(coverage.getHashes());
}
else if (hash.length() < Constants.defaultPrecision) {
coverage = GeoHash.coverBoundingBox(env.getMaxY(), env.getMinX(), env.getMinY(), env.getMaxX(), Constants.defaultPrecision);
geohashes = new ArrayList<>(coverage.getHashes());
}
else {
if (hash.length() > Constants.spaceLength) geohashes.add(hash.substring(0, hash.length() - 1));
else if (hash.length() == Constants.spaceLength) geohashes.add(hash);
else {
coverage = GeoHash.coverBoundingBox(env.getMaxY(), env.getMinX(), env.getMinY(), env.getMaxX(), hash.length() + 1);
geohashes = new ArrayList<>(coverage.getHashes());
}
}
System.out.println(geohashes.size());
List<ReferencedEnvelope> envelopes = calcuateCuts(geohashes);
ArrayList<GridCoverage2D> tiles = new ArrayList<>();
for (ReferencedEnvelope envelope : envelopes) {
tiles.addAll(getTiles(envelope, grid));
}
System.out.println(tiles.size());
return tiles;
}
它有点管用，因为我们得到了一些 tiff，但我们也得到了其中一些无法使用 qgis 等工具可视化的问题。 merge 方法给了我们一个 tiff 图像，但是它又没有在 qgis.
中可视化如果您需要更多信息，请告诉我。 谢谢！
编辑
这是原始且正确的 tif
的 gdalinfo Driver: GTiff/GeoTIFF
Files: C:\Users\Consulente\Desktop\appoggio_rilascio\originale.tif
C:\Users\Consulente\Desktop\appoggio_rilascio\originale.tif.aux.xml
Size is 266, 269
Coordinate System is:
PROJCRS["WGS 84 / UTM zone 33N",
BASEGEOGCRS["WGS 84",
DATUM["World Geodetic System 1984",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4326]],
CONVERSION["UTM zone 33N",
METHOD["Transverse Mercator",
ID["EPSG",9807]],
PARAMETER["Latitude of natural origin",0,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8801]],
PARAMETER["Longitude of natural origin",15,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Scale factor at natural origin",0.9996,
SCALEUNIT["unity",1],
ID["EPSG",8805]],
PARAMETER["False easting",500000,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",0,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["(E)",east,
ORDER[1],
LENGTHUNIT["metre",1]],
AXIS["(N)",north,
ORDER[2],
LENGTHUNIT["metre",1]],
USAGE[
SCOPE["unknown"],
AREA["World - N hemisphere - 12┬░E to 18┬░E - by country"],
BBOX[0,12,84,18]],
ID["EPSG",32633]]
Data axis to CRS axis mapping: 1,2
Origin = (293436.000000000000000,4647354.000000000000000)
Pixel Size = (3.000000000000000,-3.000000000000000)
Metadata:
AREA_OR_POINT=Area
Image Structure Metadata:
COMPRESSION=LZW
INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left ( 293436.000, 4647354.000) ( 12d30'28.04"E, 41d57' 4.05"N)
Lower Left ( 293436.000, 4646547.000) ( 12d30'29.05"E, 41d56'37.91"N)
Upper Right ( 294234.000, 4647354.000) ( 12d31' 2.66"E, 41d57' 4.80"N)
Lower Right ( 294234.000, 4646547.000) ( 12d31' 3.68"E, 41d56'38.66"N)
Center ( 293835.000, 4646950.500) ( 12d30'45.86"E, 41d56'51.36"N)
Band 1 Block=266x3 Type=UInt16, ColorInterp=Red
Min=4230.000 Max=18752.000
Minimum=4230.000, Maximum=18752.000, Mean=6493.946, StdDev=1353.641
NoData Value=0
Metadata:
STATISTICS_MAXIMUM=18752
STATISTICS_MEAN=6493.9456647057
STATISTICS_MINIMUM=4230
STATISTICS_STDDEV=1353.6413052408
STATISTICS_VALID_PERCENT=94.55
Band 2 Block=266x3 Type=UInt16, ColorInterp=Green
Min=3979.000 Max=17611.000
Minimum=3979.000, Maximum=17611.000, Mean=5971.019, StdDev=1344.279
NoData Value=0
Metadata:
STATISTICS_MAXIMUM=17611
STATISTICS_MEAN=5971.0194223549
STATISTICS_MINIMUM=3979
STATISTICS_STDDEV=1344.2792863221
STATISTICS_VALID_PERCENT=94.55
Band 3 Block=266x3 Type=UInt16, ColorInterp=Blue
Min=2619.000 Max=16659.000
Minimum=2619.000, Maximum=16659.000, Mean=4982.758, StdDev=1525.882
NoData Value=0
Metadata:
STATISTICS_MAXIMUM=16659
STATISTICS_MEAN=4982.7577379017
STATISTICS_MINIMUM=2619
STATISTICS_STDDEV=1525.8817457971
STATISTICS_VALID_PERCENT=94.55
Band 4 Block=266x3 Type=UInt16, ColorInterp=Undefined
NoData Value=0
这是来自不正确合并 tiff 的 gdalinfo（应该与原始 tiff 相同）
Driver: GTiff/GeoTIFF
Files: C:\Users\Consulente\Desktop\appoggio_rilascio\merged_tif\merged.tiff
Size is 266, 269
Coordinate System is:
GEOGCRS["WGS 84",
DATUM["World Geodetic System 1984",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433]],
CS[ellipsoidal,2],
AXIS["geodetic latitude (Lat)",north,
ORDER[1],
ANGLEUNIT["degree",0.0174532925199433]],
AXIS["geodetic longitude (Lon)",east,
ORDER[2],
ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4326]]
Data axis to CRS axis mapping: 2,1
GeoTransform =
41.95133415197179, 0, -2.77698308382438e-05
12.50778763370388, 3.722213354286106e-05, 0
Metadata:
AREA_OR_POINT=Area
TIFFTAG_RESOLUTIONUNIT=1 (unitless)
TIFFTAG_XRESOLUTION=1
TIFFTAG_YRESOLUTION=1
Image Structure Metadata:
INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left ( 41.9513342, 12.5077876) ( 41d57' 4.80"E, 12d30'28.04"N)
Lower Left ( 41.9438641, 12.5077876) ( 41d56'37.91"E, 12d30'28.04"N)
Upper Right ( 41.9513342, 12.5176887) ( 41d57' 4.80"E, 12d31' 3.68"N)
Lower Right ( 41.9438641, 12.5176887) ( 41d56'37.91"E, 12d31' 3.68"N)
Center ( 41.9475991, 12.5127382) ( 41d56'51.36"E, 12d30'45.86"N)
Band 1 Block=266x8 Type=UInt16, ColorInterp=Red
NoData Value=0
Band 2 Block=266x8 Type=UInt16, ColorInterp=Green
NoData Value=0
Band 3 Block=266x8 Type=UInt16, ColorInterp=Blue
NoData Value=0
Band 4 Block=266x8 Type=UInt16, ColorInterp=Alpha
NoData Value=0
通过查看 GdalInfo 输出，您已将 CRS 更改为 EPSG:4326 并翻转了 lat/lon 值 - 所以在您写出您没有的图块的方式中发生了一些事情发布代码。
如果我尝试这样做，我会查看 Image Tiling tutorial，关键的内部循环是这样的：
// iterate over our tile counts
for (int i = 0; i < htc; i++) {
for (int j = 0; j < vtc; j++) {
System.out.println("Processing tile at indices i: " + i + " and j: " + j);
// create the envelope of the tile
Envelope envelope = getTileEnvelope(coverageMinX, coverageMinY, geographicTileWidth, geographicTileHeight,
targetCRS, i, j);
GridCoverage2D finalCoverage = cropCoverage(gridCoverage, envelope);
if (this.getTileScale() != null) {
finalCoverage = scaleCoverage(finalCoverage);
}
// use the AbstractGridFormat's writer to write out the tile
fileExtension = "png";
File tileFile = new File(tileDirectory, i + "_" + j + "." + fileExtension);
final WorldImageWriter wiWriter = new WorldImageWriter(tileFile);
// writing parameters for png
final Format oFormat = wiWriter.getFormat();
((AbstractGridFormat) oFormat).getWriter(tileFile).write(finalCoverage, null);
}
}
和 cropCoverage 使用 CoverageProcessor class:
private GridCoverage2D cropCoverage(GridCoverage2D gridCoverage, Envelope envelope) {
CoverageProcessor processor = CoverageProcessor.getInstance();
// An example of manually creating the operation and parameters we want
final ParameterValueGroup param = processor.getOperation("CoverageCrop").getParameters();
param.parameter("Source").setValue(gridCoverage);
param.parameter("Envelope").setValue(envelope);
return (GridCoverage2D) processor.doOperation(param);
}