Lisez les données de la base de données de nuages de points de la préfecture de Shizuoka avec Java et générez une photographie aérienne et une élévation PNG.
■ Aperçu
Récemment, j'étudie le traitement de l'information géographique avec la moitié de mon travail et la moitié de mes loisirs. Lorsque je cherchais des données de groupe de points ouverts, je suis arrivé à Shizuoka Prefecture Point Cloud DB. Dans Shizuoka Prefecture Point Cloud DB, les données sont publiées dans le jeu de données Las, alors lisez Las en Java et image des données de photographie aérienne et d'élévation (altitude PNG). ) A été généré.
■ Traitement PNG d'altitude
Altitude PNG utilisé dans la tuile d'altitude du National Land Research Institute a été utilisé pour numériser l'altitude. Le papier en altitude PNG est ici. → "[PNG Altitude Tile-Consideration and Implementation of Altitude File Format Adapted for Web Use- (Yoshiharu Nishioka / Juri Nagatsu, Information Geology Vol. 26, No. 4, 2015)](https: //www.jstage.jst. go.jp/article/geoinformatics/26/4/26_155/_pdf) "
La classe ElevationPngUtil a été créée en référence au site Web de l'Institut de géographie et à l'article.
import java.awt.Color;
import java.awt.geom.AffineTransform;
import java.awt.geom.Point2D;
import java.awt.image.BufferedImage;
import java.io.BufferedWriter;
import java.io.File;
import java.io.FileWriter;
import java.io.IOException;
/**
* Classe utilitaire Altitude PNG
*/
public class ElevationPngUtil {
private static final int P8=256;
private static final int P16=65536;
private static final int P23=8388608;
private static final int P24=16777216;
private static final double U=0.01;
private ElevationPngUtil(){}
/ * Altitude PNG Valeur NA * /
public static int NA=P23;
/ * Convertit les valeurs RVB en Int * /
private static int rgb2Int(int[] c){
return new Color(c[0],c[1],c[2]).getRGB();
}
/**
* Convertir la valeur d'élévation en valeur de pixel PNG d'élévation
* @param z Valeur d'altitude (m)
* @return Altitude PNG valeur de pixel
*/
public static int getRGB(double z){
if(Double.isNaN(z))return P23;
return rgb2Int(getRGBColor(z));
}
/ * Convertit la valeur d'élévation en RVB: int [] * /
private static int[] getRGBColor(double z){
if(z<=0)return new int[]{128,0,0};
int i=(int)Math.round(z*100);
int r=i >> 16;
int g=i-(r << 16) >> 8;
int b=i-((r << 16)+(g << 8));
return new int[]{r,g,b};
}
/**
* Convertir la valeur de pixel PNG d'élévation en valeur d'élévation
* @param intColor Altitude PNG valeur de pixel
* @return Valeur d'altitude (m)
*/
public static double getZ(int intColor){
Color c=new Color(intColor);
int r=c.getRed();
int g=c.getGreen();
int b=c.getBlue();
int x=r*P16+g*P8+b;
if(x<P23){
return U*(double)x;
}else if(x>P23){
return U*(double)(x-P24);
}else{
return Double.NaN;
}
}
public static double[] getMinMaxZ(BufferedImage png){
double min=Double.MAX_VALUE;
double max=-Double.MAX_VALUE;
for(int i=0;i<png.getWidth();i++){
for(int j=0;j<png.getHeight();j++){
double h=getZ(png.getRGB(i, j));
if(Double.isNaN(h))continue;
min=Math.min(h, min);
max=Math.max(h, max);
}
}
return new double[]{min,max};
}
/**
* Convertir le PNG d'élévation en un tableau d'élévations (double [] [])
* @param png
* @return
*/
public static double[][] imgToDouble(BufferedImage png){
double[][] ret=new double[png.getWidth()][png.getHeight()];
for(int i=0;i<ret.length;i++){
for(int j=0;j<ret[i].length;j++){
ret[i][j]=getZ(png.getRGB(i, j));
}
}
return ret;
}
/**
* Convertir un tableau d'élévations (double [] []) en élévation PNG
* @param z
* @return
*/
public static BufferedImage doubleToImg(double[][] z){
BufferedImage ret=new BufferedImage(z.length,z[0].length,BufferedImage.TYPE_INT_RGB);
for(int i=0;i<z.length;i++){
for(int j=0;j<z[i].length;j++){
ret.setRGB(i, j, getRGB(z[i][j]));
}
}
return ret;
}
/**
* Soustraction d'élévation PNG
* @ param b1 Altitude PNG1
* @ param b2 Altitude PNG2
* @return Altitude PNG (Altitude PNG1-Altitude PNG2)
* @throws IllegalArgumentException
*/
public static BufferedImage subZImage(BufferedImage b1,BufferedImage b2)throws IllegalArgumentException{
int w=b1.getWidth();
int h=b1.getHeight();
if(w==b2.getWidth()&&h==b2.getHeight()){
BufferedImage ret=new BufferedImage(w,h,BufferedImage.TYPE_INT_RGB);
for(int i=0;i<w;i++){
for(int j=0;j<h;j++){
double v1=getZ(b1.getRGB(i, j));
double v2=getZ(b2.getRGB(i, j));
if(Double.isNaN(v1)||Double.isNaN(v2)){
ret.setRGB(i, j, NA);
}
double v3=v1-v2;
if(v3<0){
ret.setRGB(i, j, NA);
}else{
ret.setRGB(i, j, getRGB(v3));
}
}
}
return ret;
}else{
throw new IllegalArgumentException("Image size is different");
}
}
/**
* Sortie des données PNG d'élévation en texte
* @param png Altitude PNG
* @param af Conversion Affin
* @param out Fichier de sortie
* @throws IOException
*/
public static void output(BufferedImage png,AffineTransform af,File out) throws IOException{
Point2D dst=new Point2D.Double();
BufferedWriter bw=new BufferedWriter(new FileWriter(out));
int n=1;
int w=png.getWidth();
int h=png.getHeight();
for(int i=0;i<w;i++){
for(int j=0;j<h;j++){
dst=af.transform(new Point2D.Double(i,j), dst);
double v=getZ(png.getRGB(i, j));
if(Double.isNaN(h))continue;
bw.write(Integer.toString(n++)+","+dst.getX()+","+dst.getY()+","+Double.toString(v)+"\n");
}
bw.flush();
}
bw.close();
}
}
Chargement du jeu de données Las
J'ai utilisé "laszip4j" pour charger l'ensemble de données Las. La licence pour "laszip4j" est "GNU Lesser General Public License v2.1".
J'ai créé une classe LASFileReader qui lit et combine plusieurs fichiers Las et génère des photographies aériennes et des DSM à partir d'un groupe de points. De plus, ces données de groupe de points semblent être le système de coordonnées du 8e système à angle droit plan.
import java.awt.Color;
import java.awt.geom.AffineTransform;
import java.awt.geom.Point2D;
import java.awt.geom.Rectangle2D;
import java.awt.image.BufferedImage;
import java.io.File;
import java.io.IOException;
import java.util.ArrayList;
import java.util.Calendar;
import java.util.Date;
import java.util.List;
import javax.imageio.ImageIO;
import com.github.mreutegg.laszip4j.LASHeader;
import com.github.mreutegg.laszip4j.LASPoint;
import com.github.mreutegg.laszip4j.LASReader;
import net.termat.geo.ElevationPngUtil;
public class LASFileReader {
private LASReader reader;
private Rectangle2D bounds;
private AffineTransform af;
private BufferedImage img;
private BufferedImage dem;
private double[] minmaxZ;
private Date date;
private double xScalefactor;
private double xOffset;
private double yScalefactor;
private double yOffset;
private double zScalefactor;
private double zOffset;
private float basef=65536f;
/**
*
* @param f Las File array
* @param mPerPixel 1 longueur de pixel (m)
* @throws Exception
*/
public static void outputLasDataUnion(File[] f,double mPerPixel)throws Exception{
List<LASFileReader> list=new ArrayList<LASFileReader>();
Rectangle2D rect=null;
for(int i=0;i<f.length;i++){
if(i==0){
LASFileReader la=new LASFileReader(f[i],mPerPixel);
rect=la.getBounds();
list.add(la);
}else{
LASFileReader la=new LASFileReader(f[i],mPerPixel);
rect=rect.createUnion(la.getBounds());
list.add(la);
}
}
AffineTransform af=new AffineTransform(new double[]{mPerPixel,0,0,-mPerPixel,rect.getX(),rect.getY()+rect.getHeight()});
double width=Math.ceil(rect.getWidth())/mPerPixel;
double height=Math.ceil(rect.getHeight())/mPerPixel;
BufferedImage img=new BufferedImage((int)width,(int)height,BufferedImage.TYPE_INT_RGB);
BufferedImage dem=new BufferedImage((int)width,(int)height,BufferedImage.TYPE_INT_RGB);
for(int i=0;i<dem.getWidth();i++){
for(int j=0;j<dem.getHeight();j++){
dem.setRGB(i, j, ElevationPngUtil.NA);
}
}
for(LASFileReader al : list){
al.createImage(img,dem,af.createInverse());
}
String path=f[0].getAbsolutePath();
File pf=new File(path.replace(".las", "_photp.jpg "));
File df=new File(path.replace(".las", "_org.png "));
ImageIO.write(img, "jpg", pf);
ImageIO.write(dem, "png", df);
}
/**
*
* @param f fichier Las
* @param mPerPixel 1 longueur de pixel (m)
* @throws Exception
*/
public LASFileReader(File f,double mPerPixel)throws IOException{
reader = new LASReader(f);
readHeader();
af=new AffineTransform(new double[]{mPerPixel,0,0,-mPerPixel,bounds.getX(),bounds.getY()+bounds.getHeight()});
double width=Math.ceil(bounds.getWidth())/mPerPixel;
double height=Math.ceil(bounds.getHeight())/mPerPixel;
img=new BufferedImage((int)width,(int)height,BufferedImage.TYPE_INT_RGB);
dem=new BufferedImage((int)width,(int)height,BufferedImage.TYPE_INT_RGB);
for(int i=0;i<(int)width;i++){
for(int j=0;j<(int)height;j++){
dem.setRGB(i, j, ElevationPngUtil.NA);
}
}
}
/**
* Obtenir la zone de données Las
* @return
*/
public Rectangle2D getBounds() {
return bounds;
}
/*
* Lire l'en-tête
*/
private void readHeader(){
LASHeader h=reader.getHeader();
bounds=new Rectangle2D.Double(h.getMinX(),h.getMinY(),h.getMaxX()-h.getMinX(),h.getMaxY()-h.getMinY());
minmaxZ=new double[]{h.getMinZ(),h.getMaxZ()};
Calendar cal=Calendar.getInstance();
cal.set(Calendar.YEAR, (int)h.getFileCreationYear());
cal.set(Calendar.DAY_OF_YEAR, (int)h.getFileCreationDayOfYear());
date=cal.getTime();
xOffset=h.getXOffset();
xScalefactor=h.getXScaleFactor();
yOffset=h.getYOffset();
yScalefactor=h.getYScaleFactor();
zOffset=h.getZOffset();
zScalefactor=h.getZScaleFactor();
}
/**
* Génération d'images
* @param img photographie aérienne
* @param dem Altitude Image PNG
* @param à la conversion Affin
*/
public void createImage(BufferedImage img,BufferedImage dem,AffineTransform at){
for (LASPoint p : reader.getPoints()) {
double x=xScalefactor*(double)p.getX()+xOffset;
double y=yScalefactor*(double)p.getY()+yOffset;
double z=zScalefactor*(double)p.getZ()+zOffset;
Point2D px=at.transform(new Point2D.Double(x,y), new Point2D.Double());
float r=((float)p.getRed())/basef;
float g=((float)p.getGreen())/basef;
float b=((float)p.getBlue())/basef;
int col=new Color(r,g,b).getRGB();
img.setRGB((int)px.getX(), (int)px.getY(), col);
dem.setRGB((int)px.getX(), (int)px.getY(), ElevationPngUtil.getRGB(z));
}
}
/**
* Génération d'images
*/
public void createImage(){
AffineTransform at=null;
try{
at=af.createInverse();
}catch(Exception e){e.printStackTrace();}
for (LASPoint p : reader.getPoints()) {
double x=xScalefactor*(double)p.getX()+xOffset;
double y=yScalefactor*(double)p.getY()+yOffset;
double z=zScalefactor*(double)p.getZ()+zOffset;
Point2D px=at.transform(new Point2D.Double(x,y), new Point2D.Double());
float r=((float)p.getRed())/basef;
float g=((float)p.getGreen())/basef;
float b=((float)p.getBlue())/basef;
int col=new Color(r,g,b).getRGB();
img.setRGB((int)px.getX(), (int)px.getY(), col);
dem.setRGB((int)px.getX(), (int)px.getY(), ElevationPngUtil.getRGB(z));
}
}
}
■ Résultats
Utilisez la classe LASFileReader pour télécharger des données de groupe de points (cette fois en utilisant les données de 30XXX00010025-1.las à 30XXX00010025-6.las) depuis Shizuoka Prefecture Point Cloud DB Le traitement générera la photographie aérienne et l'image PNG d'élévation suivantes.
■ Avenir
S'il existe des données publiques telles que Shizuoka Prefecture Point Cloud DB, je pense qu'elles peuvent être utilisées pour diverses initiatives même au niveau des citoyens.
Cette fois, les données de groupe de points de la superficie forestière ont été utilisées dans "[Création de données sur la hauteur des arbres à l'aide de données de levés aériens au laser](https: //www.gsi.go.)" Introduit par le National Land Research Institute. C'est parce que je voulais essayer l'analyse de la canopée de la zone forestière en référence à "jp / chirijoho / chirijoho40069.html)".
Recommended Posts