package com.moon.server.service.data; import com.moon.server.entity.all.StaticData; import com.moon.server.entity.data.MetaEntity; import com.moon.server.helper.StringHelper; import org.apache.commons.logging.Log; import org.apache.commons.logging.LogFactory; import org.gdal.gdal.ColorTable; import org.gdal.gdal.Dataset; import org.gdal.gdal.gdal; import org.gdal.ogr.Geometry; import org.gdal.ogr.ogr; import org.gdal.osr.SpatialReference; import org.springframework.stereotype.Service; import java.io.File; /** * 栅格服务 * @author WWW * @date 2023-08-25 */ @Service public class RasterService { private SpatialReference sr4326; private SpatialReference sr4490; private SpatialReference sr104903; private final static Log log = LogFactory.getLog(RasterService.class); /** * 初始化空间引用 */ private void initSpatialReference() { if (null == sr4326) { sr4326 = new SpatialReference(); sr4326.ImportFromEPSG(StaticData.I4326); sr4490 = new SpatialReference(); sr4490.ImportFromEPSG(StaticData.I4490); sr104903 = new SpatialReference(StaticData.MOON_2000_WKT); } } /** * 读取栅格信息 */ public void readRasterInfo(MetaEntity me, String file) { Dataset ds = null; try { File f = new File(file); if (!f.exists() || f.isDirectory()) { return; } ds = gdal.Open(file, 0); if (null == ds || ds.getRasterCount() == 0) { return; } // EPSG编码 String epsg = getEpsg(ds); // 波段数 String bands = "" + ds.getRasterCount(); //vd.Meta.band_type = Enum.GetName(typeof(DataType), ds.GetRasterBand(1).DataType); // 数据类型 ColorTable colorTable = ds.GetRasterBand(1).GetRasterColorTable(); // 数据颜色表 String ct = null == colorTable ? null : colorTable.toString(); // 高程基准 String hDatum = null; double[] transform = new double[6]; ds.GetGeoTransform(transform); // 分辨率 String resolution = String.format("%f,%f", transform[1], transform[5]); if (!StaticData.EPSGS.contains(epsg)) { return; } Geometry minPoint = getMinPoint(ds); Geometry maxPoint = getMaxPoint(ds); double xmin = minPoint.GetX(0); double ymin = minPoint.GetY(0); double xmax = maxPoint.GetX(0); double ymax = maxPoint.GetY(0); // 四至范围 String geom = String.format("ST_GeomFromText('POLYGON ((%f %f,%f %f,%f %f,%f %f,%f %f))')", xmin, ymax, xmax, ymax, xmax, ymin, xmin, ymin, xmin, ymax); } catch (Exception ex) { log.error(ex.getMessage(), ex); } finally { if (null != ds) { ds.delete(); } } } /** * 获取EPSG编码 */ private String getEpsg(Dataset ds) { SpatialReference sr = ds.GetSpatialRef(); if (sr != null) { // 坐标系统 String coorSys = sr.GetName(); if (StaticData.MOON200.equals(coorSys)) { return "ESRI:" + StaticData.I104903; } else { String code = sr.GetAuthorityCode(null); return StringHelper.isEmpty(code) ? null : "EPSG:" + code; } } return null; } /** * 获取Dataset的最小点 */ private Geometry getMinPoint(Dataset ds) { double[] transform = new double[6]; ds.GetGeoTransform(transform); String epsg = ds.GetSpatialRef().GetAuthorityCode(null); double xMin = transform[0]; double yMin = transform[3] - ds.getRasterYSize() * transform[1]; Geometry point = new Geometry(ogr.wkbPoint); point.AddPoint(xMin, yMin, 0); return transform(ds, point, epsg); } /** * 获取Dataset的最大点 */ private Geometry getMaxPoint(Dataset ds) { double[] transform = new double[6]; ds.GetGeoTransform(transform); String epsg = ds.GetSpatialRef().GetAuthorityCode(null); double xMax = transform[0] + (ds.getRasterXSize() * transform[1]); double yMax = transform[3]; Geometry point = new Geometry(ogr.wkbPoint); point.AddPoint(xMax, yMax, 0); return transform(ds, point, epsg); } /** * 坐标转换 */ private Geometry transform(Dataset ds, Geometry point, String epsg) { this.initSpatialReference(); if (StringHelper.isEmpty(epsg)) { point.AssignSpatialReference(sr104903); return point; } if (String.valueOf(StaticData.I4326).equals(epsg)) { point.AssignSpatialReference(sr4326); return point; } if (String.valueOf(StaticData.I4490).equals(epsg)) { point.AssignSpatialReference(sr4490); return point; } point.AssignSpatialReference(ds.GetSpatialRef()); if (ds.GetSpatialRef().GetName().contains(StaticData.CGCS2000)) { point.TransformTo(sr4490); } else { point.TransformTo(sr4326); } return point; } }