Suite

Traiter le raster : le créer dans le CRS projeté et l'exporter dans le CRS géographique

Traiter le raster : le créer dans le CRS projeté et l'exporter dans le CRS géographique


Lors de la création d'un raster à partir d'un polygone (le SCR du polygone est un système de coordonnées géographiques), le raster est créé dans le système de référence de coordonnées projetées du bloc de données, qui est un système projeté. Je dois le faire car je dois spécifier la taille de la cellule du raster en mètres.

Je dois également livrer ce raster (précédemment converti en polygone) en WGS84 (GCS); J'envisage de le faire: Changer le système de référence de coordonnées du bloc de données en WGS84, et, sur le raster de résultat, faire Données> Exporter les données et, pour l'option "Utiliser le même système de coordonnées que:", choisissez "le bloc de données ”.

Est-ce une façon précise de traiter cela? Je veux dire, n'y a-t-il pas des modifications dans les zones ou les coordonnées de la couche de forme livrée ?


EDIT : Lorsque vous reprojetez votre raster il y aura un rééchantillonnage, donc les nouvelles données seront dégradées. Il est donc préférable de convertir directement à partir d'une classe d'entités avec le CRS en sortie qui est le plus adapté à votre besoin, en d'autres termes vous devez reprojeter votre classe d'entités qui est actuellement dans WGS 84 en un CRS projeté si vous voulez que vos données finales soient dans WGS 84. La conversion de la classe d'entités sera exactement réversible si vous ne modifiez pas le datum, vous devez juste être conscient que seuls les sommets sont projetés. par conséquent, si vous avez de longs segments sans sommets, il est recommandé de "[densifier][1]" avant votre projet.

La taille de vos pixels affectera la précision de la conversion de votre entité en raster, ainsi que la précision du rééchantillonnage. Cela n'affectera cependant pas la précision de la surface de votre polygone liée à la projection (si votre projection n'est pas à surface égale, la sur ou sous-estimation moyenne de la surface des polygones ne dépendra pas de la taille du pixel)

Considérons maintenant le choix d'un SIR.

Avec un système de coordonnées géographiques (lat/long), la taille de la cellule sera exprimée en degré. Les cellules d'un raster de la même largeur et de la même hauteur partout, donc l'étendue réelle (au sol) des cellules de la grille changera en fonction de votre latitude. Si vous voulez connaître la taille en mètre d'une cellule dans un système lat/long, vous pouvez avoir une bonne approximation en sachant qu'un degré vaut environ 111,3 km à l'équateur, mais que la largeur de votre pixel doit être corrigée par le cosinus de la latitude. Par conséquent, si vous avez besoin d'un raster final en Lat/long, oubliez une taille de cellule exacte en mètre et utilisez une approximation pratique basée sur votre emplacement (par exemple, le minimum ou la moyenne d'un degré de Lat et d'un degré de Long). Par exemple, à une latitude de 60°N, un carré d'un mètre sur un mesure environ 8,9847e-6 degré dans la direction Y et 1,797e-5 degré dans la direction X.

Avec un système de coordonnées projetées, la taille de la cellule est exprimée dans un système cartésien (avec une unité en mètre, en pieds…) et cette valeur est constante partout sur la grille. Bien sûr, vous êtes dans un système de coordonnées projetées, il y a donc une petite distorsion due au passage d'une sphère à un plan, mais vous ne devriez pas trop avoir peur de cette distorsion car vous pouvez généralement trouver des projections locales qui réduisent la distorsion en dessous de la plage d'erreur de votre outil/méthode de mesure. Par conséquent, si vous n'avez pas besoin d'avoir un raster en latitude/longitude, utilisez un système de coordonnées projetées (local) dans lequel vous projetez votre classe d'entités, puis convertissez la classe d'entités projetée en raster.

Remarque : Si vous n'effectuez pas votre analyse en géométrie sphérique, je vous recommande de travailler dans un SIR projeté.


QGIS prend en charge environ 2 700 CRS connus. Les définitions de chaque SCR sont stockées dans une base de données SQLite installée avec QGIS. Normalement, vous n'avez pas besoin de manipuler directement la base de données. En fait, cela peut entraîner l'échec du support de projection. Les CRS personnalisés sont stockés dans une base de données utilisateur. Voir section Système de référence de coordonnées personnalisé pour plus d'informations sur la gestion de vos systèmes de référence de coordonnées personnalisés.

Les SIR disponibles dans QGIS sont basés sur ceux définis par l'European Petroleum Search Group (EPSG) et l'Institut Géographique National de France (IGNF) et sont largement abstraits des tables de référence spatiales utilisées dans GDAL. Les identifiants EPSG sont présents dans la base de données et peuvent être utilisés pour spécifier un SCR dans QGIS.

Pour utiliser la projection OTF, soit vos données doivent contenir des informations sur son système de référence de coordonnées, soit vous devrez définir un CRS global, de couche ou à l'échelle du projet. Pour les couches PostGIS, QGIS utilise l'identifiant de référence spatiale qui a été spécifié lors de la création de la couche. Pour les données prises en charge par OGR, QGIS s'appuie sur la présence d'un moyen reconnu de spécifier le CRS. Dans le cas des fichiers de formes, cela signifie un fichier contenant la spécification de texte bien connue (WKT) du CRS. Ce fichier de projection a le même nom de base que le fichier de formes et un .prj extension. Par exemple, un fichier de formes nommé alaska.shp aurait un fichier de projection correspondant nommé alaska.prj .

Chaque fois que vous sélectionnez un nouveau SCR, les unités de calque seront automatiquement modifiées dans le Général onglet du Propriétés du projet dialogue sous le Projet (Gnome, OS X) ou Paramètres (KDE, Windows).


Projection cartographique en détail¶

Une méthode traditionnelle pour représenter la forme de la terre est l'utilisation de globes. Il y a cependant un problème avec cette approche. Bien que les globes conservent la majorité de la forme de la terre et illustrent la configuration spatiale des entités de la taille d'un continent, ils sont très difficiles à transporter dans une poche. Ils ne sont également pratiques à utiliser qu'à des échelles extrêmement petites (par exemple 1:100 million).

La plupart des données cartographiques thématiques couramment utilisées dans les applications SIG sont à une échelle considérablement plus grande. Les jeux de données SIG typiques ont des échelles de 1:250 000 ou plus, selon le niveau de détail. Un globe de cette taille serait difficile et coûteux à produire et encore plus difficile à transporter. En conséquence, les cartographes ont développé un ensemble de techniques appelées projections cartographiques conçu pour montrer, avec une précision raisonnable, la terre sphérique en deux dimensions.

Vu de près, la terre semble relativement plate. Cependant, vu de l'espace, nous pouvons voir que la terre est relativement sphérique. Les cartes, comme nous le verrons dans le prochain sujet sur la production de cartes, sont des représentations de la réalité. Ils sont conçus pour représenter non seulement les caractéristiques, mais aussi leur forme et leur disposition spatiale. Chaque projection cartographique a avantages et désavantages. La meilleure projection pour une carte dépend de la échelle de la carte et sur les fins pour lesquelles elle sera utilisée. Par exemple, une projection peut avoir des distorsions inacceptables si elle est utilisée pour cartographier l'ensemble du continent africain, mais peut être un excellent choix pour un carte à grande échelle (détaillée) de votre pays. Les propriétés d'une projection cartographique peuvent également influencer certaines des caractéristiques de conception de la carte. Certaines projections sont bonnes pour de petites zones, certaines sont bonnes pour cartographier des zones avec une grande étendue est-ouest, et certaines sont meilleures pour cartographier des zones avec une grande étendue nord-sud.


Ensembles de données

GeoTiff est un format d'image typique pour les données raster. Si vous n'avez pas de fichier GeoTiff à portée de main, vous pouvez en télécharger un ici. Les exemples suivants sont basés sur ce fichier. Nous pouvons lire dans un fichier GeoTiff dans un base de données, la structure de données principale de rasterio, en utilisant le code suivant :

data est appelé un objet de jeu de données, qui contient des informations utiles, notamment :

EPSG:32631 est un système de référence local utilisé pour certaines parties de l'Europe. http://epsg.io fournit des informations sur la plupart des systèmes de référence de coordonnées disponibles. Les coordonnées Crs sont généralement étiquetées x et y dans le cas de EPSG:4326 (qui est le même que WGS84) les coordonnées sont fournies en degrés dans la définition commune de longitude et de latitude.

Les nombres fournis utilisent les unités comme crs de l'ensemble de données dans ce cas spécifique (EPSG:32631), les nombres définissent la boîte englobante de la carte en unités de mètres par rapport à une origine définie par les crs. D'autres définitions de crs peuvent utiliser d'autres unités, par exemple les degrés.

Cette transformation, implémentée en tant qu'objet Affine, définit comment un changement de 1 pixel dans l'une ou l'autre direction (ligne ou colonne) se traduit par des changements de coordonnées crs en utilisant 6 paramètres qui sont (dans cet ordre et tous dans les mêmes unités que les crs de la transformation) : le changement de x en fonction du changement de colonne de pixels (+10m pour +1 pixel), le changement de x en fonction du changement de rangée de pixels (0 dans ce cas), l'origine de la coordonnée x (590520.0m) , le changement de y en fonction du changement de colonne de pixels (0 dans ce cas), le changement de y en fonction du changement de rangée de pixels (-10m pour +1 pixel), et l'origine de la coordonnée y (5790630.0 Je suis là). Vous pouvez considérer que la transformation consiste en deux parties : une « matrice de rotation » et l'origine des coordonnées crs du jeu de données (voir le croquis en haut de cet article).

(Dans ce cas, 3 bandes différentes sont disponibles, chaque bande représente une carte en niveaux de gris pour une région de longueur d'onde spécifique. Notez que rasterio commence à compter les bandes à 1, contrairement à Python, pour une raison quelconque.)

Nous pouvons lire une bande spécifique de l'ensemble de données dans un tableau en utilisant

où i est l'un des indices fournis par data.indexes . img n'est en réalité qu'un tableau numpy, notez que tous les tableaux de bandes disponibles ont la même forme (et donc également l'orientation et la mise à l'échelle) car ils partagent les mêmes crs, transformation et limites.

Traçons notre exemple d'image. Nous pouvons utiliser les méthodes rasterio.plot Adjust_band , qui normalise les valeurs du tableau dans une plage de zéro à un, et show , qui crée et affiche les axes matplotlib. Notez que nous fournissons à imgdata l'ordre inverse de l'index (ou du canal) puisque l'ensemble de données contient les canaux BGR, mais nous voulons tracer RVB :


Volet Exporter raster

Le volet Exporter le raster vous permet d'exporter l'intégralité du jeu de données raster, de la mosaïque, du service d'imagerie ou la partie de l'affichage.

Fournit des options pour spécifier le jeu de données raster en sortie , le système de coordonnées , les transformations géographiques , la géométrie de découpage , le maintien de l'étendue de découpage , la taille de cellule , la taille de raster , le type de pixel , la valeur de pixel d'échelle , la valeur NoData , les paramètres de rendu , le format de sortie , le type de compression et la qualité de compression .

Vous permet de configurer le raster de capture , la taille des tuiles , le rééchantillonnage , le type de source , les paramètres des pyramides et les paramètres des statistiques.

Vous pouvez également ouvrir ce volet en cliquant sur le bouton Exporter le raster dans l'onglet Données.

Lors du stockage d'un jeu de données raster dans une géodatabase, aucune extension de fichier ne doit être ajoutée au nom du jeu de données raster. Lors du stockage du jeu de données raster dans un format de fichier, vous devez spécifier l'extension appropriée au fichier.

La boîte de dialogue Référence spatiale est contextuelle et répertorie un système de coordonnées xy ou z selon que vous sélectionnez l'option Système de coordonnées XY actuel ou Z actuel.

  1. Sélectionnez la transformation géographique appropriée lorsque vos données sont transformées entre différents systèmes de coordonnées. L'application n'utilisera que les transformations appropriées à la projection, toutes les autres seront ignorées.

Cette option exportera le jeu de données raster à l'aide des spécifications de référencement spatial du jeu de données raster.

L'étendue de l'affichage actuel sera utilisée.

Par exemple, si vous effectuez un zoom avant sur votre zone d'étude particulière, vous pouvez utiliser cette option pour traiter les entités qui se trouvent dans l'étendue d'affichage actuelle.

Vous entrez les coordonnées du rectangle de délimitation minimum : saisissez l'étendue pour Gauche , Droite , Haut et Bas .

Toutes les couches sont répertoriées et vous pouvez en choisir une à utiliser comme étendue.

Comme avec l'option Étendue d'affichage actuelle, l'étendue de la couche est lue et stockée.

Utilisez le bouton Parcourir pour rechercher l'emplacement du dossier de la classe d'entités que vous souhaitez utiliser pour la géométrie de découpage . Une fois qu'une entité en entrée est fournie, une case à cocher apparaît avec l'option Utiliser les entités en entrée pour la géométrie de découpage, avec des options de découpage pour délimiter l'intérieur ou l'extérieur .

Conserve l'alignement des cellules en tant que raster en entrée et ajuste l'étendue en sortie en conséquence.

Ajuste le nombre de colonnes et de lignes et rééchantillonne les pixels pour qu'ils correspondent exactement à l'étendue de découpage spécifiée.

L'option Mettre à l'échelle la valeur de pixel apparaît lorsqu'un type de pixel différent est choisi, qui peut être utilisé pour mettre à l'échelle votre type de pixel d'une profondeur de bits à une autre.

Si le type de pixel est rétrogradé (abaissé), les valeurs raster en dehors de la plage valide pour cette profondeur de pixel seront tronquées et perdues. Par exemple, lorsque la sortie est un type de pixel différent de l'entrée (comme 16 bits à 8 bits), vous pouvez choisir de mettre les valeurs à l'échelle pour s'adapter à la nouvelle plage sinon, les valeurs qui ne rentrent pas dans le nouveau pixel la plage sera supprimée.

En cas de mise à l'échelle, telle que 8 bits à 16 bits, le minimum et le maximum des valeurs 8 bits seront mis à l'échelle au minimum et au maximum dans la plage 16 bits. En cas de réduction, telle que 16 bits à 8 bits, le minimum et le maximum des valeurs 16 bits seront réduits au minimum et au maximum dans la plage 8 bits. Pour en savoir plus sur la capacité de profondeur de bits pour les formats d'exportation pris en charge, voir Liste des capteurs pris en charge.

Les valeurs de pixel resteront les mêmes et ne seront pas mises à l'échelle. Toutes les valeurs qui ne rentrent pas dans la plage de valeurs seront rejetées. C'est la valeur par défaut.

Les valeurs de pixel seront mises à l'échelle au nouveau type de pixel. Lorsque vous mettez à l'échelle la profondeur de vos pixels, votre raster affichera la même chose, mais les valeurs seront mises à l'échelle à la nouvelle profondeur de bits qui a été spécifiée.

Ceci est recommandé si vous exportez vers un jeu de données raster basé sur un fichier et que l'écrêtage graphique est choisi. Lorsqu'un graphique est utilisé pour découper vos données, des pixels NoData existeront très probablement dans la sortie. La spécification de la valeur NoData vous permet de contrôler la profondeur de pixel et la valeur qui stockera NoData, cependant, la valeur des pixels NoData n'est pas stockée pour les données exportées vers une géodatabase ou un format CRF.

Cochez la case Forcer RVB pour exporter le raster en sortie en tant que jeu de données raster RVB à trois canaux avec le moteur de rendu actuel.

De plus, lors de l'exportation vers des formats tels que TIFF, JP2, PNG et MRF qui prennent en charge les bandes alpha, vous pouvez utiliser cette option pour exporter les données en tant que jeu de données raster à quatre canaux avec une bande alpha pour préserver les paramètres de transparence des données d'origine .

La case à cocher Utiliser la palette de couleurs est activée uniquement si votre raster source contient une palette de couleurs ou si le rendu de couche raster actuel est un rendu à valeur unique.

Cochez la case Utiliser le rendu pour exporter le jeu de données raster avec les statistiques et les options actuelles du rendu. Lorsque vous ouvrez le jeu de données raster exporté dans ArcGIS Pro , les règles de rendu par défaut sont appliquées. Pour conserver le même rendu que lorsque vous avez exporté les données, définissez le type d'étirement sur Aucun , car il est déjà étiré.

  • La case à cocher Utiliser le rendu doit être cochée pour utiliser Force RVB ou Utiliser la palette de couleurs .
  • Les options Force RGB et Use Colormap ne peuvent pas être utilisées ensemble.
  • La combinaison de Use Renderer et Force RGB exporte la sortie en tant que raster RVB avec trois ou quatre bandes (bande alpha le cas échéant) avec l'affichage du rendu de couche raster actuel lors de l'exportation.
  • La combinaison de Use Renderer et Use Colormap exporte le raster avec une palette de couleurs dans les cas suivants :
    • Si le raster source contient une palette de couleurs.
    • Pour les jeux de données qui ne contiennent pas de palette de couleurs, mais peuvent être visualisés dans un moteur de rendu à valeur unique.
    1. Choisissez le type de compression si votre format de sortie le permet.
    2. Choisissez la qualité de compression si votre format de sortie est JP2 ou JPG.

    Lorsque les rasters sont stockés sous forme de blocs de données, les jeux de données raster sont stockés dans un type de données appelé BLOB (binary large object). L'option de taille de tuile vous permet de contrôler le nombre de pixels stockés dans chaque BLOB et, par conséquent, vous permet de contrôler la taille de chaque BLOB. Il est spécifié comme le nombre de pixels dans X (largeur de la tuile) et Y (hauteur de la tuile).

    Le rééchantillonnage est le processus d'interpolation des valeurs de pixels tout en transformant votre jeu de données raster. Ceci est utilisé lorsque l'entrée et la sortie ne s'alignent pas exactement, lorsque la taille des pixels change, lorsque les données sont décalées, ou une combinaison de ceux-ci.

    Effectue une affectation de voisin le plus proche et est la plus rapide des méthodes d'interpolation. Il est principalement utilisé pour les données discrètes, telles qu'une classification de l'utilisation des terres, car il ne modifiera pas les valeurs des pixels. L'erreur spatiale maximale sera d'un demi-pixel.

    Effectue une interpolation bilinéaire et détermine la nouvelle valeur d'un pixel en fonction d'une moyenne de distance pondérée des quatre centres de pixel d'entrée les plus proches. Il est utile pour les données continues et entraînera un certain lissage des données.

    Effectue une convolution cubique et détermine la nouvelle valeur d'un pixel en fonction de l'ajustement d'une courbe lisse passant par les 16 centres de pixels d'entrée les plus proches. Il convient aux données continues, bien qu'il puisse en résulter que le raster en sortie contienne des valeurs en dehors de la plage du raster en entrée. Il est géométriquement moins déformé que le raster obtenu en exécutant l'algorithme de rééchantillonnage du voisin le plus proche. L'inconvénient de l'option Cubic est qu'elle nécessite plus de temps de traitement. I Si le temps de traitement est un problème, utilisez plutôt Bilinéaire.

    Il n'y a pas de type de données spécifié.

    Le raster est un type de données d'altitude.

    Le raster est un type de données thématique qui a des valeurs discrètes, telles que la couverture terrestre.

    La trame a été traitée en couleur et ne doit pas être étirée en contraste.

    Le raster contient des informations scientifiques et sera affiché avec la rampe de couleur bleu à rouge par défaut.

    Le raster est un raster à deux canaux qui contient une composante U et une composante V des données de champ vectoriel.

    Le raster est un raster à deux canaux qui contient l'amplitude et la direction des données de champ vectoriel.

    1. Spécifiez le nombre de niveaux de pyramide . Vous pouvez spécifier le nombre de niveaux à créer, ou vous pouvez laisser la valeur vide pour créer tous les niveaux.
    2. Cochez la case Ignorer le premier pour ignorer le premier niveau de pyramide de votre raster.
    3. Spécifiez la technique de rééchantillonnage de la pyramide utilisée pour construire vos pyramides : Voisin le plus proche , Bilinéaire ou Cubique .
    4. Choisissez le type de compression de la pyramide à utiliser lors de la création des pyramides raster.
    • Par défaut : le système détectera un type de compression approprié. Si les données source sont compressées à l'aide d'une compression par ondelettes, elle construira des pyramides avec le type de compression JPEG sinon, LZ77 sera utilisé.
    • Aucun : aucune compression n'est utilisée lors de la création de pyramides.
    • LZ77 : l'algorithme de compression LZ77 sera utilisé pour construire les pyramides. LZ77 peut être utilisé pour tout type de données.
    • JPEG : l'algorithme de compression JPEG sera utilisé pour créer des pyramides. Seules les données conformes à la spécification de compression JPEG peuvent utiliser ce type de compression. Si JPEG est choisi, vous pouvez alors définir la qualité de compression .
    • JPEG YCbCr : une compression avec perte utilisant les composants de l'espace colorimétrique luma (Y) et chroma (Cb et Cr).
    1. Choisissez si vous souhaitez sauter des pixels entre les échantillons. Les paramètres X Skip Factor et Y Skip Factor représentent respectivement le nombre de pixels horizontaux et verticaux entre les échantillons. La valeur doit être supérieure à zéro et inférieure ou égale au nombre de colonnes ou de lignes du jeu de données raster.
    2. Le paramètre Statistiques ignorer les valeurs vous permet d'ignorer une ou plusieurs valeurs qui ne participeront pas au calcul des statistiques, comme une valeur d'arrière-plan. Les valeurs multiples sont séparées par des points-virgules.

    Une fois l'exportation terminée, le jeu de données raster exporté est ajouté à la carte.


    Natural Earth a un jeu de données Admin 0 - Pays. Télécharger les pays

    L'Ordnance Survey du Royaume-Uni fournit des données ouvertes à télécharger. Téléchargez le produit raster MiniScale pour la Grande-Bretagne et extrayez-le dans un dossier sur votre ordinateur.

    Vous devrez entrer vos informations personnelles pour pouvoir télécharger l'ensemble de données Ordnance Survey.

    Pour plus de commodité, vous pouvez télécharger directement une copie de l'ensemble de données à partir du lien ci-dessous :

    minisc_gb.zip (Contient uniquement les fichiers requis pour ce tutoriel)


    Natural Earth a un jeu de données Admin 0 - Pays. Télécharger les pays

    L'Ordnance Survey du Royaume-Uni fournit des données ouvertes à télécharger. Téléchargez le produit raster MiniScale pour la Grande-Bretagne et extrayez-le dans un dossier sur votre ordinateur.

    Vous devrez entrer vos informations personnelles pour pouvoir télécharger l'ensemble de données Ordnance Survey.

    Pour plus de commodité, vous pouvez télécharger directement une copie de l'ensemble de données à partir du lien ci-dessous :

    minisc_gb.zip (Contient uniquement les fichiers requis pour ce tutoriel)


    Récupération de données raster par emplacement géographique à l'aide de Landsat et PostGIS

    Le projet sur lequel je travaille nécessite que je récupère des données raster Landsat à des emplacements géographiques spécifiques (long/lat). Après avoir parcouru quelques tutoriels et expérimenté avec GDAL, PostGIS et QGIS, j'ai importé avec succès une image GeoTIFF Landsat dans une table raster PostGIS et j'ai accédé aux valeurs par emplacement géographique à partir de cette table. Cependant, il y avait quelques problèmes dans le résultat:

    • Je ne comprends pas le système de coordonnées utilisé par QGIS dans son interface, car ils se situent dans les centaines de milliers
    • Le raster chargé dans QGIS au large des côtes espagnoles, plutôt qu'au-dessus du Maine, aux États-Unis, comme il était censé le faire.

    Voici quelques informations sur ma démarche. Je suis assez nouveau sur les SIG en général, donc je suis presque certain qu'il y a une erreur flagrante à trouver ici :

    • Télécharger le fichier Landsat 8 GeoTIFF de l'USGS GloVis
    • Renommez l'image du groupe 5 en quelque chose de plus convivial pour commander un ninja.
    • Créez une base de données postgres pour les tables raster et exécutez CREATE EXTENSION postgis

    Exécutez gdalinfo LSSampleB5.TIF , en imprimant la sortie suivante :

    Pilote : GTiff/GeoTIFF Fichiers : LSSampleB5Test2.TIF La taille est 7871, 7971 Le système de coordonnées est : PROJCS["WGS 84 / UTM zone 19N", GEOGCS["WGS 84", DATUM["WGS_1984", SPHEROID["WGS 84", 6378137,298.257223563, AUTORITE["EPSG","7030"]], AUTORITE["EPSG","6326"]], PRIMEM["Greenwich",0, AUTORITE["EPSG","8901"]], UNITÉ[ "degré",0.0174532925199433, AUTORITE["EPSG","9122"]], AUTORITE["EPSG","4326"]], PROJECTION["Transverse_Mercator"], PARAMETRE["latitude_of_origin",0], PARAMETRE["central_meridian ",-69], PARAMETER["scale_factor",0.9996], PARAMETER["false_easting",500000], PARAMETER["false_northing",0], UNIT["metre",1, AUTHORITY["EPSG","9001" ]], AXIS["Easting",EAST], AXIS["Northing",NORTH], AUTHORITY["EPSG","32619"]] Origine = (318285,0000000000000000000,216715,000000000000000) Taille de pixel = (30.0000000000000000,-30.0000000000000000) Métadonnées : AREA_OR_POINT=Métadonnées de la structure de l'image de point : INTERLEAVE=BAND Coordonnées du coin : en haut à gauche ( 318285.000, 5216715.000) ( 71d23'37.53"W, 47d 4'44.12"N) en bas à gauche ( 318 285.000, 4977585.000) ( 71d18' 9.77"W, 44d55'42.53"N) En haut à droite ( 554415.000, 5216715.000) ( 68d16'58.41"W, 47d 6' 6.11"N) En bas à droite ( 554415.000, 4977585.000) ( 68d18'36.69" W, 44d56'58.62"N) Centre ( 436350.000, 5097150.000) ( 69d49'20.56"W, 46d 1'29.87"N) Bande 1 Bloc=7871x1 Type=UInt16, ColorInterp=Gray

    J'ai interprété cette sortie comme le format EPSG 4326 (ce qui peut être mon crime), j'ai donc exécuté la commande suivante pour importer le GeoTIFF en tant que raster PostGIS :

    raster2pgsql -s 4326 -I LSSampleB5.TIF -F -t 50x50 -d | psql -U postgres test raster

    Cela a importé avec succès une nouvelle table. J'ai ensuite utilisé QGIS pour avoir une intuition visuelle de ce qui se passait.

    Sous Base de données -> Gestionnaire de base de données -> PostGIS -> rastertest -> public, j'ai ajouté mon lssampeb5 au canevas.

    J'ai créé une nouvelle connexion XYZ dans QGIS pour ajouter des images hybrides satillite de Google à titre de référence. L'url que j'ai utilisé était https://mt1.google.com/vt/lyrs=y&x=&y=&z= avec un zoom min et max de 0 et 19 respectivement.

    C'est ici que j'ai pris note du fait que la couche d'échantillon a atterri au large des côtes espagnoles sur la carte Google Hybrid.

    Je me suis assuré que les deux couches étaient sur la projection EPSG 4326, aucun changement.

    Pas trop découragé pour passer à autre chose, j'ai essayé une requête de base de données pour obtenir une seule valeur de pixel. Étant donné que mes données d'échantillon ont atterri près de l'Espagne, j'ai utilisé QGIS pour échantillonner une paire de coordonnées valide près de là pour la requête. La requête était :

    SELECT rid, ST_Value(rast, 1, ST_SetSRID(ST_Point(448956,5041439), 4326)) as b5 FROM lssampeb5 WHERE ST_Intersects(rast, ST_SetSRID(ST_Point(448956,5041439), 4326)::geometry, 1)

    Cela a renvoyé un ID de ligne valide et une ST_VALUE de 5776. Essayer des coordonnées en dehors de la plage affichée par QGIS n'a entraîné aucune entrée renvoyée, ce qui n'est pas inattendu.

    Donc, tout d'abord, je ne sais pas ce que QGIS utilise pour son système de coordonnées. Ce n'est certainement pas la longitude et la latitude sous une forme brute, mais d'après ma compréhension, EPSG 4326 est censé être une projection géographique.

    Deuxièmement, je ne sais pas pourquoi QGIS égare la scène Landsat au mauvais endroit, ou à quel endroit dans le processus la scène n'a pas été transformée correctement.


    Exemple reproductible

    SF points objet spatial

    En bref : le CRS est principalement stocké au format WKT. L'ancien proj4string est disponible sur demande et ne laisse pas tomber la partie datum/towgs84

    SP points objet spatial

    En bref : le CRS est principalement stocké au format chaîne proj4 et supprime la partie datum et towgs84 contrairement à sf ). la nouvelle notation WKT est stockée en commentaire dans l'objet CRS mais est différente de sf

    Raster

    Dans raster , l'ancienne notation proj4 et la nouvelle notation WKT semblent être stockées.

    Cet ensemble de données ne contient pas la partie +towgs84 dans la chaîne proj4. Mais lorsque vous lisez une trame avec une partie + towgs84 dans la chaîne proj4, elle semble être supprimée.
    Exemple non reproductible :

    Je devrais probablement aussi explorer ce qui se passe lorsque nous utilisons le package stars au lieu de raster mais cette question est déjà assez longue (et j'ai beaucoup de code construit sur le package raster)


    Étape 5- Recadrage

    Maintenant que nous avons une grande carte d'occurrence d'eau fusionnée dans le même CRS que notre image Landsat 8, il est temps de la recadrer en utilisant les étendues de notre zone d'intérêt. Donc, la première chose à faire maintenant est d'obtenir l'étendue de l'image Landsat. Nous avons vu dans Étape 2 que nous pouvons utiliser la méthode des limites pour obtenir l'étendue d'un ensemble de données. Cependant, la méthode de recadrage de rasterio reçoit un GeoJSON en entrée pour effectuer le recadrage.

    Nous pouvons extraire le GeoJSON de notre image landsat en utilisant le rasterio.features.shapes fonction, qui nous renverra un GeoJSON des formes qu'il identifie dans un ensemble de données donné. Comme nous voulons l'image entière, nous pouvons passer juste un tableau vide avec la même forme que l'image Landsat. L'astuce ici est que les fonctions de formes renvoient un "générateur", qui nous permet d'itérer à travers les différents résultats. Nous n'entrerons pas dans les détails du rendement et des générateurs ici, mais comme nous n'attendons qu'une seule sortie (le polygone avec l'image entière), nous allons utiliser le python intégré Suivant fonction, comme ceci :

    Comme avant, le recadrer La fonction nécessite un ensemble de données comme argument, nous allons donc créer le merged_ds au préalable. Ensuite, nous recadrerons dans les limites que nous venons de récupérer et afficherons les résultats.