Partie 1 Si vous cliquez sur le bouton du plug-in J'ai créé un plug-in qui peut convertir les coordonnées de la couche multi-polygone EPSG: 4301 en EPSG: 4612.
Cependant, il est inutile uniquement pour les couches multi-polygones. Puisqu'il ne prend pas en charge les couches dans le système de coordonnées orthogonales plan Cette fois, je voudrais renforcer ce domaine et le rendre un peu plus pratique.
・ Tous les types de géométrie peuvent être convertis ・ Peut être converti indépendamment de la latitude / longitude et de l'angle du plan -Une couche de post-conversion qui hérite de la valeur d'attribut peut être créée. ~~ ・ Vous pouvez sélectionner plusieurs couches et les convertir à la fois ~~ ・ Ne gérez pas les exceptions
Étant donné que la conversion de plusieurs couches à la fois est une question liée à l'interface utilisateur, j'ai décidé de les séparer. Faisons de notre mieux.
・ Conversion avant (système de géographie japonais → système de géographie mondiale) Ce que vous pouvez faire (terminé) </ font> ・ Peut être converti indépendamment de la latitude / longitude et de l'angle du plan ・ Tous les types de géométrie peuvent être convertis -Une couche de post-conversion qui hérite de la valeur d'attribut peut être créée. ・ Interpole la valeur de correction du fichier de paramètres pour améliorer la précision. * Difficulté: moyenne ・ Plusieurs couches peuvent être sélectionnées et converties à la fois ・ La conversion inverse (système d'enquête mondiale → système d'enquête japonais) est également possible * Difficulté: élevée
Si vous n'avez pas créé la source dans la partie 1, supprimez-la du référentiel. https://github.com/ozo360/tky2jgd/releases/tag/ddbc2ba
Je ne dis pas que cela correspond à tout. Le type de géométrie de la couche vectorielle que QGIS peut gérer est QgsWkbTypes :: geometryType Un autre type est QgsWkbTypes :: Type.
La cible de ce plug-in est limitée aux données 2D avec des points, des lignes et des polygones. Cela devrait couvrir la plupart des données. Je ne sais pas s'il peut être converti correctement avec les données Curve et Triangle, je vais donc le tester après l'avoir fait. Si cela ne fonctionne pas, modifiez-le pour l'exclure.
Pour obtenir la couche source de conversion wkbType et l'intégrer dans l'URI de la nouvelle couche, procédez comme suit. Le traitement des exceptions n'est pas effectué, mais une branche est incluse pour limiter la cible.
layer = self.iface.activeLayer()
wkbType = layer.wkbType()
#polygone de ligne cible
if not QgsWkbTypes.geometryType(wkbType) in [
QgsWkbTypes.PointGeometry,
QgsWkbTypes.LineGeometry,
QgsWkbTypes.PolygonGeometry]:
#Non applicable geometryType
self.iface.messageBar().pushMessage("tky2jgd <b>{}</b>".format('activeLayer is out of scope geometryType'), Qgis.Warning)
#Cibler les couches 2D
if QgsWkbTypes.hasZ(wkbType) or QgsWkbTypes.hasM(wkbType):
#Non applicable wkbType
self.iface.messageBar().pushMessage("tky2jgd <b>{}</b>".format('activeLayer is out of scope wkbType'), Qgis.Warning)
self.wkbType = QgsWkbTypes.displayString(wkbType)
layerUri = self.wkbType +"?crs=postgis:4612"
newlayer = QgsVectorLayer(layerUri, "exchange", "memory")
Cela m'a pris du temps car je ne savais pas comment obtenir la chaîne d'énumération (displayString), mais j'ai réussi à créer une couche qui correspond au type de géométrie de la couche cible.
La cible étant le système géodésique japonais, il est nécessaire de pouvoir effectuer les conversions suivantes.
système | SRID de la source de conversion | SRID de destination de conversion |
---|---|---|
4301 | 4612 | |
I | 30161 | 2443 |
II | 30162 | 2444 |
III | 30163 | 2445 |
IV | 30164 | 2446 |
・ | ・ | ・ |
・ | ・ | ・ |
・ | ・ | ・ |
XVIII | 30178 | 2460 |
XIX | 30179 | 2461 |
En outre, le fichier de paramètres affiche la valeur de correction dans EPSG: 4301. L'angle droit du plan doit être converti en latitude et longitude avant de déplacer la valeur de correction. Par conséquent, la conversion vers le SRID de destination de conversion sera à partir d'EPSG: 4301 pour tout SRID. Générez une classe de conversion de coordonnées pour chacun.
layer = self.iface.activeLayer()
srid = layer.crs().postgisSrid()
if srid == 4301:
self.isLonlat = True
self.toSrid = 4612
elif (srid > 30160 and srid < 30180):
self.isLonlat = False
self.toSrid = 2442 + srid - 30160
else:
#Sans objet SRID
self.iface.messageBar().pushMessage("tky2jgd <b>{}</b>".format('activeLayer is not Tokyo Datum'), Qgis.Warning)
#Générer une classe de transformation de coordonnées
crs4301 = QgsCoordinateReferenceSystem(4301, QgsCoordinateReferenceSystem.EpsgCrsId)
if not self.isLonlat:
crsFrom = QgsCoordinateReferenceSystem(fromSrid, QgsCoordinateReferenceSystem.EpsgCrsId)
self.trans4301 = QgsCoordinateTransform(crsFrom, crs4301, QgsProject.instance())
crsTo = QgsCoordinateReferenceSystem(self.toSrid, QgsCoordinateReferenceSystem.EpsgCrsId)
self.transTo = QgsCoordinateTransform(crs4301, crsTo, QgsProject.instance())
Tant que vous avez une classe de conversion de coordonnées, s'il s'agit d'un plan orthogonal au moment du traitement de conversion, il vous suffit d'ajouter un prétraitement pour le convertir en latitude et longitude.
Vous devez ajouter une colonne à la nouvelle couche pour hériter des valeurs d'attribut. Copiez les colonnes du calque d'origine vers le nouveau calque.
#Quand newLayer est un nouveau calque
layer = self.iface.activeLayer()
newLayer.dataProvider().addAttributes(layer.fields())
Lors de la création d'une fonction, la valeur d'attribut d'origine est héritée.
inFeat = QgsFeature()
feats = layer.getFeatures()
while feats.nextFeature( inFeat ):
attributes = inFeat.attributes()
geometry = QgsGeometry( inFeat.geometry() )
#En fait, les coordonnées sont converties ici
feature = QgsFeature(fields)
feature.setAttributes(attributes)
feature.setGeometry(geometry)
newLayer.addFeature(feature)
Maintenant que chaque processus est terminé, nous allons l'intégrer dans le plug-in.
###################################
#Cliquez sur le bouton de démarrage du processus
###################################
def buttonClicked(self):
self.layer = self.iface.activeLayer()
if not self.isScope():
return
self.wkbType = QgsWkbTypes.displayString(self.layer.wkbType())
#Générer un dictionnaire de paramètres de conversion de coordonnées
self.loadPar()
#Générer une classe de transformation de coordonnées
self.setCrsTrans()
#Processus de conversion
self.execTrans()
QMessageBox.information(self.dlg, 'tky2jgd', 'finished')
###################################
#Lire le fichier de paramètres de conversion
###################################
def loadPar(self):
self.par = {}
parfile = 'TKY2JGD.par'
with open(os.path.join(os.path.abspath(os.path.dirname(__file__)), parfile)) as f:
#Passer deux lignes dans l'en-tête
skip = 2
for i in range(skip):
next(f)
#Lire ligne par ligne et stocker dans la liste
while True:
line = f.readline()
if not line:
# EOF
break
#La valeur est en secondes, alors divisez-la ici
self.par[int(line[0:8])] = (float(line[8:18]) / 3600, float(line[18:28]) / 3600)
###################################
#Déterminer s'il s'agit d'un calque cible de conversion
#(Au fait, conservez les indicateurs srid et latitude / longitude avant et après la conversion)
###################################
def isScope(self):
try:
if not isinstance(self.layer, QgsVectorLayer):
raise ValueError("activeLayer is not QgsVectorLayer")
self.srid = self.layer.crs().postgisSrid()
if self.srid == 4301:
self.isLonlat = True
self.toSrid = 4612
elif (self.srid > 30160 and self.srid < 30180):
self.isLonlat = False
self.toSrid = 2442 + self.srid - 30160
else:
#Sans objet SRID
raise ValueError("activeLayer is not Tokyo Datum")
wkbType = self.layer.wkbType()
#polygone de ligne cible
if not QgsWkbTypes.geometryType(wkbType) in [
QgsWkbTypes.PointGeometry,
QgsWkbTypes.LineGeometry,
QgsWkbTypes.PolygonGeometry]:
#Non applicable geometryType
raise ValueError("activeLayer is out of scope geometryType")
#Cibler les couches 2D
if QgsWkbTypes.hasZ(wkbType) or QgsWkbTypes.hasM(wkbType):
#Non applicable wkbType
raise ValueError("activeLayer is out of scope wkbType")
except ValueError as e:
#N'est pas applicable
self.iface.messageBar().pushMessage("tky2jgd <b>{}</b>".format(e), Qgis.Warning)
return False
else:
#Cible
return True
###################################
#Coordonner le processus de conversion
###################################
def execTrans(self):
#Générer une couche après la conversion
afterLayer = self.createAfterLayer()
#Commencer l'édition
afterLayer.startEditing()
fields = afterLayer.fields()
inFeat = QgsFeature()
feats = self.layer.getFeatures()
while feats.nextFeature( inFeat ):
attributes = inFeat.attributes()
beforeGeom = QgsGeometry( inFeat.geometry() )
if not self.isLonlat:
beforeGeom.transform(self.trans4301)
afterGeom = self.moveCorrection(beforeGeom)
if not self.isLonlat:
afterGeom.transform(self.transTo)
feature = QgsFeature(fields)
feature.setAttributes(attributes)
feature.setGeometry(afterGeom)
afterLayer.addFeature(feature)
#Fin de l'édition
afterLayer.commitChanges()
QgsProject.instance().addMapLayer(afterLayer)
###################################
#Créer une couche de mémoire convertie
###################################
def createAfterLayer(self):
layerUri = self.wkbType +"?crs=postgis:" + str(self.toSrid)
afterLayer = QgsVectorLayer(layerUri, "exchange", "memory")
afterLayer.dataProvider().addAttributes(self.layer.fields())
return afterLayer
###################################
#Déplacement du haut d'une entité
###################################
def moveCorrection(self, geom):
for i, v in enumerate(geom.vertices()):
meshcode = self.Coordinate2MeshCode(v.y(), v.x())
correction = self.par[meshcode]
geom.moveVertex(v.x() + correction[1], v.y() + correction[0], i)
return geom
###################################
#Générer une classe de transformation de coordonnées
###################################
def setCrsTrans(self):
crs4301 = QgsCoordinateReferenceSystem(4301, QgsCoordinateReferenceSystem.EpsgCrsId)
if not self.isLonlat:
crsFrom = QgsCoordinateReferenceSystem(self.srid, QgsCoordinateReferenceSystem.EpsgCrsId)
self.trans4301 = QgsCoordinateTransform(crsFrom, crs4301, QgsProject.instance())
crsTo = QgsCoordinateReferenceSystem(self.toSrid, QgsCoordinateReferenceSystem.EpsgCrsId)
self.transTo = QgsCoordinateTransform(crs4301, crsTo, QgsProject.instance())
Je n'ai pas commencé correctement avec la conception des classes, donc c'était un peu compliqué. Quant à cet espace, nous allons le découper en classes à la fin et l'organiser, donc nous donnerons la priorité au déménagement en premier. C'est tellement ennuyeux que la fonction isScope me donne envie de mourir.
Il est difficile de créer des échantillons de données sobrement.
Vérifions que les données d'échantillon créées par nous-mêmes peuvent être converties même à angle droit par rapport au plan.
Le SRID de la couche convertie est correctement 2458. © OpenStreetMap contributors
couche | X | Y |
---|---|---|
Avant la conversion | -21165 | -177520 |
Après la conversion | -21313 | -177066 |
Version Web TKY2JGD | -21315.4052 | -177083.9327 |
L'erreur dans la direction Y est grande, mais si le calcul est faux, cette erreur ne devrait pas suffire, alors peut-être que ça va?
Continuons à l'essayer à Tokyo. En même temps, le type de géométrie est vérifié.
Les couches multipoints peuvent également être converties et le SRID de cette couche est 2451. © OpenStreetMap contributors
couche | X | Y |
---|---|---|
Avant la conversion | -6984 | -34730 |
Après la conversion | -7276 | -34371 |
Version Web TKY2JGD | -7277.1503 | -34374.4408 |
Je n'ai pas de mauvais pressentiment, mais je ferme les yeux et passe à autre chose.
Maintenant, vérifiez si les attributs sont hérités Ajoutez des colonnes de caractères, d'entiers et de virgules flottantes à la couche mémoire, puis ajoutez des entités. Donnons une valeur d'attribut. Lorsque les attributs du calque converti sont affichés, ils sont correctement conservés.
Avec cela, toutes les fonctions qui faisaient l'objet de cette époque ont été implémentées.
Cliquez ici pour la source jusqu'à présent
Essayez d'enregistrer le calque converti dans un Shapefile. Bien entendu, le CRS spécifié est celui du calque converti. Lorsque je charge à nouveau le Shapefile enregistré dans QGIS, il change. Je dois me rendre à la volée près du calque d'origine, mais j'arrive à l'emplacement converti. Lorsque je vérifie les coordonnées, elles ne sont vraiment pas alignées. Une valeur similaire à la valeur obtenue en convertissant les anciennes coordonnées en nouvelles coordonnées, puis en les déplaçant vers les nouvelles coordonnées.
Vous avez peut-être fait une mauvaise façon de l'ajouter à votre projet, ou vous avez peut-être commis une erreur fondamentale avant cela. Normalement, je pense que vous devriez publier après avoir résolu le problème, mais il est normal d'avoir de tels essais et erreurs. C'est pourquoi je vais dormir, donc je ne peux pas le mettre à jour pendant un moment, mais la prochaine fois, j'étudierai ce bogue.
Après avoir déplacé avec le fichier de paramètres dans la conversion de la latitude et de la longitude, c'était une erreur de convertir davantage les anciennes coordonnées en nouvelles coordonnées avec QgsCoordinateTransform, donc je l'ai corrigé en incluant "Part 1". En conséquence, la conversion de la latitude et de la longitude est devenue normale. Il a été constaté que la quantité de mouvement du fichier de paramètres était doublée lors de la conversion de l'orthogonalité du plan. La cause de ceci n'a pas encore été trouvée.
![Licence Creative Commons](https://i.creativecommons.org/l/by/4.0 /88x31.png) Cet article est fourni sous la Creative Commons Attribution 4.0 International License .
Recommended Posts