[PYTHON] Devinons l'état de développement de la ville à partir de l'image satellite.
Aperçu
Jusqu'à présent, nous avons acquis des images d'observation du ** Copernicus hub **, qui est une plaque tournante pour les images satellites Sentinel, et avons obtenu des tendances concernant l'état de stationnement des parkings à Tokyo Disneyland et le nombre de navires dans la baie de Tokyo.
J'ai vérifié l'état d'utilisation du parking à partir d'images satellite. J'ai essayé de trouver la tendance du nombre de navires dans la baie de Tokyo à partir d'images satellites.
Cette fois-ci, comme par le passé, nous obtiendrons des images satellites de ** Google Earth Engine ** et évaluerons l'état de développement de la banlieue de Tokyo au cours des cinq dernières années à partir d'images satellites.
Veuillez consulter le site ci-dessus pour savoir comment utiliser Google Earth Engine (ci-après, GEE). ** Si vous avez un compte Google, vous pouvez l'utiliser avec Google colaboratory **, vous n'avez donc pas à créer un environnement difficile, vous pouvez donc en faire l'expérience relativement facilement **.
Il s'agit d'une image composite de l'image d'observation ** de Sentinel-1 le 12 septembre 2020 ** et de l'image d'observation ** du 21 septembre 2015 ** il y a environ 5 ans. L'image passée est rouge et l'image récente est bleue.
De là, des parties bleues relativement grandes peuvent être vues à Odaiba et à Tokyo Disneyland. On estime que ces lieux ont été ** développés au cours des 5 dernières années **.
Je présenterai la méthode de création en incluant l'explication de cette image composite. La sortie finale crée un fichier kmz pour visualiser les images satellites ci-dessus sur Google Earth. Ce fichier se trouve à ici, donc si vous voulez l'essayer sur votre Google Earth, veuillez le télécharger. De plus, le code introduit dans l'article est placé sur Github, donc si vous voulez l'essayer vous-même, veuillez l'utiliser. ** Google colaboratory est utilisé **, vous pouvez donc l'utiliser quel que soit le PC tant que vous disposez d'un environnement réseau **.
1.Tout d'abord
Cette fois, à partir des images satellites, nous utiliserons les images prises par le ** satellite SAR ** appelé capteur d'ondes radio. L'image capturée du satellite SAR est ** difficile à structurer **, mais il peut être évalué qualitativement s'il existe de nombreuses structures, et elle est stable car elle n'est pas affectée par les nuages. peut faire. Cependant, il y a une différence de résolution même avec la même image SAR, et j'aimerais utiliser une image SAR haute résolution si possible, mais cette fois aussi ** l'image d'observation du satellite Sentinel-1 qui est libre et régulièrement observée ** Je vais l'utiliser. La densité des structures par les images SAR peut être déduite dans une certaine mesure des images. Je n'ai pas de cœur de dessin, je vais donc me référer aux articles connexes.
Kihon de données satellitaires - Compréhension, type, fréquence, résolution, cas d'utilisation
~
Le satellite SAR irradie la Terre avec des ondes radio et des images de l'intensité du signal réfléchi. À ce moment, ** les ondes radio sont émises à partir d'une direction oblique par rapport à la surface du sol **, donc si la surface du sol est plane comme indiqué sur la figure, beaucoup de choses sont dans la direction incidente des ondes radio, c'est-à-dire le satellite qui reçoit les ondes radio émises. Il réfléchit dans la direction opposée. Par conséquent, ** Si la surface au sol cible est plate, le signal reçu est faible et l'image sera sombre. La surface de l'eau des lacs et de la mer est extrêmement plate, elle semble donc plus sombre **. D'autre part, lorsqu'il y a une structure telle qu'un bâtiment, un bâtiment ou une voiture, l'objet est inégal, donc plus d'ondes radio sont réfléchies dans la même direction que la direction incidente, de sorte que l'image semble lumineuse. Par exemple, l'image d'observation du satellite radio de Tokyo Disneyland est la suivante.
Vous pouvez voir les installations de Disneyland. De plus, comme vous le savez probablement si vous êtes allé à Disneyland, vous pouvez voir qu'il y a un parking autour, qui est large et plat, donc il fait sombre. En regardant la page d'accueil de Disneyland, il est prévu de développer un projet d'expansion à grande échelle de Tokyo Disney Sea à l'endroit qui était un parking en 2015. ， Il s'est avéré que l'attraction est en cours de construction ici en ce moment. En outre, Odaiba et la région du Golfe sont en cours de développement en fonction d'événements tels que les Jeux olympiques de Tokyo, et vous pourrez peut-être voir le modèle à partir d'images satellite. Les satellites artificiels utilisés cette fois sont Sentinel-1 et 2 développés et exploités par l'Agence spatiale européenne. Pour plus de détails sur chacun, veuillez consulter les articles suivants.
Comment obtenir gratuitement les dernières images satellite.
L'image observée de l'évaluation d'analyse est Sentinel-1, qui est un satellite SAR.
Sentinel-1 est composé de deux satellites, ** l'un a 12 jours de retour ** (c'est-à-dire qu'il traverse la même zone dans les mêmes conditions tous les 12 jours), mais il est composé de deux. Par conséquent, nous observons de nombreuses zones ** une fois tous les 6 jours **. Cette ** observation régulière ** est une caractéristique majeure du satellite Sentinel et un avantage côté utilisateur.
[Sentinel-1 Observation Scenario]
(https://sentinel.esa.int/web/sentinel/missions/sentinel-1/observation-scenario)
De plus, l'heure d'observation des satellites artificiels est quasiment fixe **, pour les satellites optiques vers 10h30, et pour les satellites SAR vers 6h00 et 18h00 **.
Pouvez-vous voir des gens depuis le satellite artificiel? ~ Par satellite, résolution au sol / résumé de l'heure locale ~
Pour Sentinel-1, qui a été introduit plus tôt, l'heure d'observation est d'environ 6 heures du matin pour le satellite A et de 18 heures environ pour le satellite B. Étant donné que la direction d'irradiation des ondes radio est différente pour chaque satellite, il n'est pas approprié de comparer en utilisant les deux images satellite car la différence peut être capturée. Par conséquent, cette fois, nous n'utiliserons que l'image du satellite B pour acquérir l'image passée à comparer avec l'image satellite récente et en créer une image composite. Le but de la création d'une image composite est de rendre les modifications visuellement faciles à comprendre. En rendant l'image récente en bleu et l'image passée en rouge, la partie qui apparaît en bleu indique que le signal est devenu plus fort, et la partie qui apparaît en rouge indique que le signal est devenu plus faible. Comme la force du signal dépend de la présence ou de l'absence de bâtiments et de la densité, ** la partie bleue indique que le développement a progressé **, et ** la partie rouge indique que la structure a disparu ** (aplatie). Il est. Ensuite, nous présenterons l'acquisition d'images satellites de GEE et son évaluation.
2. Acquisition et évaluation d'images satellites
2.1 Préparation de l'environnement (construction).
Je l'ai présenté jusqu'à présent, mais comme certains d'entre vous sont nouveaux dans cet article, je vais le dupliquer. Si vous savez, veuillez passer à 2.4. Veuillez vous référer à l'article suivant pour plus de détails sur l'acquisition d'images satellite avec GEE. Analyse d'images satellites artificielles par Google Earth Engine et Google Colab-Analyse d'images satellite à partir de gratuitement (Introduction) - ** Facile à utiliser gratuitement si vous avez un compte Google **. C'est devenu un monde incroyable. Les données acquises telles que les images satellite sont enregistrées dans Google Drive. Étant donné que le code et les données sont supprimés chaque fois que vous utilisez Google Colaboratory, il est pratique de conserver les données dans votre propre Google Drive. Cependant, étant donné que la capacité de Google dirve pour une utilisation gratuite est de 20 Go, elle s'épuisera dès qu'une grande image satellite sera téléchargée, alors veillez à la supprimer le cas échéant. Maintenant, laissez-moi vous présenter le code.
import ee
import numpy as np
import matplotlib.pyplot as plt
ee.Authenticate()
ee.Initialize()
Tout d'abord, exécutez ceci pour authentifier la connexion GEE. Lorsque vous l'exécutez, le lien sera retourné. Cliquez dessus pour effectuer la procédure d'authentification, copiez le code d'accès et saisissez-le.
Ensuite, authentifiez la connexion à Google Drive. Encore une fois, le flux est le même que celui de la certification GEE.
from google.colab import drive
drive.mount('/content/drive')
Ensuite, nous effectuerons des travaux tels que la visualisation de l'image satellite acquise et l'installation des modules nécessaires à sa numérisation et à son analyse.
#Installation du package&importer
!pip install rasterio
import numpy as np
import matplotlib.pyplot as plt
import rasterio
import json
import os
import glob
import time
from datetime import datetime
from dateutil.parser import parse
Les modules fréquemment utilisés sont déjà installés dans Google Colaboratory, donc aucun travail supplémentaire n'est requis, mais cette fois, nous utiliserons Geotiff, qui est une image avec des informations cartographiques ajoutées, il est donc nécessaire pour le traitement d'image ** Rasterio Installez **.
Ensuite, installez un module appelé ** folium ** pour vérifier la zone cible définie sur la carte.
!pip install folium
import folium
Maintenant que l'environnement est prêt, prenez l'image satellite de GEE.
2.2 Définition de la zone d'intérêt (zone cible)
Pour obtenir des images satellite de GEE, vous devez entrer ** les informations de latitude / longitude ** de la zone cible qui vous intéresse. Je me souviens en quelque sorte de la latitude et de la longitude du Japon lors de la classe de géographie à l'école et dans le monde, mais quelle est la latitude et la longitude de ma maison? Même si c'est dit, je pense que vous pouvez enfin le découvrir en cherchant et en cherchant. (Bien sûr, je ne me souviens pas non plus.) Par conséquent, j'ai fait ce qui suit afin que la latitude et la longitude de la zone d'intérêt puissent être facilement étudiées et obtenues.
#Acquisition des informations polygonales de la zone d'intérêt.
from IPython.display import HTML
HTML(r'<iframe width="1000" height="580" src="https://gispolygon.herokuapp.com/" frameborder="0"></iframe>')
Lorsque vous exécutez cela, l'écran suivant s'affiche.
Après avoir élargi la zone d'intérêt, sélectionnez le polygone carré dans l'icône de gauche pour afficher le polygone de la zone d'intérêt. Après cela, cliquez sur ** Afficher les entités ** pour afficher les informations géographiques du polygone dans la fenêtre de droite. Cliquez ensuite sur ** Copier ** en bas pour copier ces informations géographiques.
Par exemple, dans le cas du parking de Tokyo Disneyland, qui est la cible de cette fois, ce sera comme suit.
Collez ensuite les informations cartographiques copiées ci-dessous et entrez-les.
A = {"type":"FeatureCollection","features":[{"properties":{"note":"","distance":"56641.31 m","drawtype":"rectangle","area":"38085.35 ha"},"type":"Feature","geometry":{"type":"Polygon","coordinates":[[[139.7144973278046,35.51642373317183],[139.7144973278046,35.6731886200871],[139.95559573173526,35.6731886200871],[139.95559573173526,35.51642373317183],[139.7144973278046,35.51642373317183]]]}}]}
Ces informations géographiques sont traitées dans le format d'entrée pour GEE et pour l'affichage dans Folium.
#Définissez n'importe quel nom de fichier à utiliser ultérieurement. Par exemple, le nom de la zone.
object_name = 'TokyoBay_RGB'
with open(str(object_name) +'_2.geojson', 'w') as f:
json.dump(A, f)
json_file = open(str(object_name) +'_2.geojson')
json_object = json.load(json_file)
#Seules les informations de latitude / longitude de la zone d'intérêt sont extraites de json.
AREA = json_object["features"][0]["geometry"]['coordinates'][0]
area = pd.DataFrame(AREA, columns=['longtitude', 'latitude'])
area_d =[[area['longtitude'].min(), area['latitude'].max()],
[area['longtitude'].max(), area['latitude'].max()],
[area['longtitude'].max(), area['latitude'].min()],
[area['longtitude'].min(), area['latitude'].min()],
[area['longtitude'].min(), area['latitude'].max()]]
AREA = area_d
Maintenant, vérifions la zone d'intérêt définie.
m = folium.Map([(AREA[0][1]+AREA[len(AREA)-2][1])/2,(AREA[0][0]+AREA[len(AREA)-3][0])/2], zoom_start=11)
folium.GeoJson(str(object_name) +'_2.geojson').add_to(m)
m
production
2.3 Acquisition d'images satellites de GEE
De nombreuses images satellites et de nombreuses informations déjà analysées sont définies dans GEE. Pour plus de détails, veuillez consulter Data Catalog. Sentinel-1 et 2 sont les suivants.
Sentinel-1 SAR GRD: C-band Synthetic Aperture Radar Ground Range Detected, log scaling Sentinel-2 MSI: MultiSpectral Instrument, Level-1C
À partir de cette page, ** les images d'observation Sentinel-1 peuvent être utilisées à partir du 3 octobre 2014 **, et ** les données Sentinel-2 peuvent être utilisées à partir du 23 juin 2015 **. Quant à l'image d'observation de Sentinel-2, une image de niveau 2A de correction atmosphérique avec le brouillard de l'atmosphère enlevé est également préparée, mais seule l'image postérieure au 28 mars 2017 lorsque l'analyse est devenue standard est ciblée. Il devient. Par conséquent, si vous souhaitez utiliser une image d'observation plus ancienne, veuillez utiliser cette image.
Maintenant, récupérez les images de Sentinel-1 et 2 de GEE et enregistrez-les dans Google Colaboratory.
Tout d'abord, préparez le format des informations géographiques à définir dans GEE.
region=ee.Geometry.Rectangle(area['longtitude'].min(),area['latitude'].min(), area['longtitude'].max(), area['latitude'].max())
Ensuite, définissez les paramètres des informations à acquérir. Cette fois, la période de l'image acquise et la destination de sauvegarde de l'image acquise sont spécifiées. Cette fois, afin d'acquérir des images satellites de jours différents, nous allons d'abord acquérir les dernières images satellites. Ici, c'est le 20 septembre 2020.
#Précisez la date d'échéance
import datetime
date = datetime.date(2020, 9, 20)
dt_now = datetime.date.today()
dt = dt_now - date
if dt.days < 6:
dt_nd = datetime.timedelta(days=12)
dt_pd = datetime.timedelta(days=0)
else:
dt_nd = datetime.timedelta(days=6)
dt_pd = datetime.timedelta(days=6)
#Période d'observation
from_date= date - dt_nd
to_date= date + dt_pd
from_date = from_date.strftime('%Y-%m-%d')
to_date = to_date.strftime('%Y-%m-%d')
#Nom du dossier à enregistrer
dir_name_s1 = 'GEE_Sentinel1_' + object_name
dir_name_s2 = 'GEE_Sentinel2_' + object_name
Ici, puisque la fréquence d'observation de Sentinel-1 est une fois tous les 12 jours, si la date définie est proche de l'heure d'analyse, l'image passée sera acquise, et s'il s'agit d'une certaine image précédente, l'image pendant 6 jours avant et après cela sera acquise. fait.
Maintenant, définissons les conditions d'image pour Sentinel-1 et 2. Cette fois, je n'utiliserai pas l'image de Sentinel-2, mais je l'obtiendrai pour référence.
def cloudMasking(image):
qa = image.select('QA60')
cloudBitMask = 1 << 10
cirrusBitMask = 1 << 11
mask = qa.bitwiseAnd(cloudBitMask).eq(0).And(qa.bitwiseAnd(cirrusBitMask).eq(0))
return image.updateMask(mask).divide(10000)
def ImageExport(image,description,folder,region,scale):
task = ee.batch.Export.image.toDrive(image=image,description=description,folder=folder,region=region,scale=scale)
task.start()
#Seules les images ascendantes ou décroissantes sont utilisées pour comparer les images observées sur le même satellite.
Sentinel1 = ee.ImageCollection('COPERNICUS/S1_GRD').filterBounds(region).filterDate(parse(from_date),parse(to_date)).filter(ee.Filter.eq('orbitProperties_pass', 'ASCENDING')).select(['VV'])
Sentinel2 = ee.ImageCollection('COPERNICUS/S2').filterBounds(region).filterDate(parse(from_date),parse(to_date)).filterMetadata('CLOUDY_PIXEL_PERCENTAGE','less_than', 10).map(cloudMasking).select(['B4','B3','B2'])
imageList_s1 = Sentinel1.toList(300)
imageList_s2 = Sentinel2.toList(300)
Ici, puisque l'image d'observation de Sentinel-1 n'utilise que l'image d'observation du satellite ** B à 18h00 **, ** 'ASCENDING' ** est sélectionné dans 'orbitProperties_pass'. Si vous le réglez sur "Décroissant", l'image d'observation sera à 6h00.
Désormais, dans les conditions d'acquisition d'images satellite ci-dessus, les images des zones d'intérêt de Sentinel-1 et Sentinel-2 seront acquises.
for i in range(imageList_s1.size().getInfo()):
image = ee.Image(imageList_s1.get(i))
ImageExport(image.reproject(crs='EPSG:4326',scale=10),image.get('system:index').getInfo(),dir_name_s1,region['coordinates'][0],10)
for i in range(imageList_s2.size().getInfo()):
image = ee.Image(imageList_s2.get(i))
ImageExport(image.reproject(crs='EPSG:4326',scale=10),image.get('system:index').getInfo(),dir_name_s2,region['coordinates'][0],10)
2.4 Affichage des images satellites et traitement de lissage par filtres
Maintenant, affichons et confirmons l'image satellite acquise en 2.3. Les images satellites sont enregistrées dans le répertoire (dossier) défini dans mon lecteur de Google Drive. Appelez-le et affichez-le. Tout d'abord, à partir de l'image d'observation de Sentinel-2.
#Visualisation par ordre chronologique
s2_path = '/content/drive/My Drive/' + dir_name_s2 + '/'
files =os.listdir(s2_path)
files.sort()
plt.figure(figsize=(25, 25))
j=0
v = len(files)//5 +1
for i in range(len(files)):
#Acquérir et visualiser des images une scène à la fois
with rasterio.open(s2_path + files[i]) as src:
arr = src.read()
j+=1#Décaler et placer la position de tracé de l'image
plt.subplot(v,5,j)
arrayImg = np.asarray(arr).transpose(1,2,0).astype(np.float32)*2 #La luminosité est corrigée en doublant la luminosité.
plt.imshow(arrayImg)
plt.title(files[i][0:8])#Obtenir la date à partir du nom du fichier
#plt.tight_layout()
** La résolution de l'image observée de Sentinel-2 est de 10 m **, mais s'il s'agit d'une grande structure, vous savez peut-être si elle a été construite. Étant donné que les images d'observation des satellites optiques nous sont familières sur Google Map, nous pouvons facilement les trouver car le but (l'emplacement) est spécifié lors de l'utilisation de ce service. Cependant, cela peut être remarqué s'il y a des informations selon lesquelles ** la structure y aurait été construite **, mais découvrez si la structure est construite quelque part dans cette plage. , Sera difficile ** car il n'y a pas d'index **. Par conséquent, cette fois, l'image d'observation de Sentinel-2 n'est pas utilisée pour l'évaluation, mais en montrant le changement en comparant l'image d'observation de Sentinel-1, ** il est possible de trouver facilement la différence **. Maintenant, affichons l'image observée de Sentinel-1.
#Visualisation par ordre chronologique
s1_path = '/content/drive/My Drive/' + dir_name_s1 + '/'
files =os.listdir(s1_path)
files.sort()
plt.figure(figsize=(20, 40))
j=0
v = len(files)//5 +1
for i in range(len(files)):
#Acquérir et visualiser des images une scène à la fois
with rasterio.open(s1_path + files[i]) as src:
arr = src.read()
#print(arr[0].shape)
j+=1#Décaler et placer la position de tracé de l'image
plt.subplot(v,5,j)
plt.imshow(arr[0], cmap='gray')
plt.title(files[i][33:41])#Obtenir la date à partir du nom du fichier
date0 = files[i][33:41]
plt.tight_layout()
À partir de là, deux images d'observation avec la même date d'observation ont été acquises. On pense que la raison pour laquelle les deux images d'observation sont séparées dans la même zone est que la taille standard de l'image satellite de Sentinel-1 est de 100 km ☓ 100 km, et la zone d'intérêt l'a traversée cette fois.
Par conséquent, les images satellites divisées en deux sont d'abord combinées en une seule.
#Lire les données
n = 0
with rasterio.open(s1_path + files[n]) as src:
arr = src.read()
print(files[n][33:41])
#Visualisation
plt.imshow(arr[0], cmap='gray')
type(arr[0])
image1 =np.zeros(arr[0].shape)
image_1 = np.nan_to_num(arr[0])
#Lire les données
n = 1
with rasterio.open(s1_path + files[n]) as src:
arr = src.read()
print(files[n][33:41])
#Visualisation
plt.imshow(arr[0], cmap='gray')
type(arr[0])
image2 =np.zeros(arr[0].shape)
image_2 = np.nan_to_num(arr[0])
Ici, la zone (partie blanche) qui n'a pas été acquise dans la zone d'intérêt est ** (Nan) ** sans données, elle est donc remplacée par ** 0 **. Ensuite, combinez ces deux images.
image_BG = image_1 + image_2
plt.imshow(image_BG, cmap='gray')
Vous avez maintenant préparé l'image satellite du satellite SAR le 12 septembre 2020. L'image du satellite SAR produit un grand nombre de bruits de points lumineux appelés ** bruit de speckle ** en raison des réflexions des structures au sol. Ceci est causé par des interférences lors de la réflexion car les ondes radio utilisées par le satellite SAR sont des ondes. ** Un traitement de lissage ** est nécessaire car ce bruit peut affecter le traitement ultérieur. Cette fois, nous avons utilisé le ** filtre Lee **, qui est un processus général de lissage d'image SAR. La fonction de filtre Lee préparée est la suivante.
from scipy.ndimage.filters import uniform_filter
from scipy.ndimage.measurements import variance
def lee_filter(img, size):
img_mean = uniform_filter(img, (size, size))
img_sqr_mean = uniform_filter(img**2, (size, size))
img_variance = img_sqr_mean - img_mean**2
overall_variance = variance(img)
img_weights = img_variance / (img_variance + overall_variance)
img_output = img_mean + img_weights * (img - img_mean)
return img_output
Maintenant, effectuons le processus de lissage. Ici, le processus de lissage à moins de 3 pixels a été effectué.
image_BG_lee = lee_filter(image_BG, 3)
plt.imshow(image_BG_lee,cmap='gray')
Vous ne pouvez pas voir l'effet du filtre avec la taille de cette image satellite. Veuillez vérifier la différence en téléchargeant l'image créée. Ensuite, comparez les images satellites passées. Ici, l'image a été prise il y a 5 ans.
#Spécifiez la date d'observation de l'image de base
date = datetime.date(2015, 9, 20)
print(date)
#Période d'observation
from_date= date - dt_nd
to_date= date + dt_pd
from_date = from_date.strftime('%Y-%m-%d')
to_date = to_date.strftime('%Y-%m-%d')
#Nom du dossier à enregistrer
dir_name_s1 = 'GEE_Sentinel1_R_' + object_name
dir_name_s2 = 'GEE_Sentinel2_R_' + object_name
for i in range(imageList_s1.size().getInfo()):
image = ee.Image(imageList_s1.get(i))
ImageExport(image.reproject(crs='EPSG:4326',scale=10),image.get('system:index').getInfo(),dir_name_s1,region['coordinates'][0],10)
Maintenant, affichons l'image acquise.
#Visualisation par ordre chronologique
s1_path = '/content/drive/My Drive/' + dir_name_s1 + '/'
files =os.listdir(s1_path)
files.sort()
plt.figure(figsize=(20, 40))
j=0
v = len(files)//5 +1
for i in range(len(files)):
#Acquérir et visualiser des images une scène à la fois
with rasterio.open(s1_path + files[i]) as src:
arr = src.read()
print(arr[0].shape)
#print(image)
j+=1#Décale la position du tracé de l'image
plt.subplot(v,5,j)
#image = arr[0].fillna(0)
plt.imshow(arr[0], cmap='gray')
plt.title(files[i][33:41])#Obtenir la date à partir du nom du fichier
date1 = files[i][33:41]
plt.tight_layout()
Cette fois, l'image différente de la précédente n'est pas divisée, et c'est une image comprenant toute la zone d'intérêt. Cette image est également lissée avec un ** filtre ** comme auparavant.
image_R_lee = lee_filter(image_BG, 3)
plt.imshow(image_R_lee,cmap='gray')
Ensuite, créez une ** image RVB ** qui combine la dernière image satellite et l'image satellite précédente. Cette fois, ** les images passées sont définies en rouge ** et ** les images récentes sont en bleu et vert **.
2.5 Créer une image composite RVB d'images satellites récentes et passées
Maintenant, créez une ** image composite RVB ** à partir de chaque image satellite.
#Composition RVB
RGB = np.dstack((image_R, np.dstack((image_BG_lee, image_BG_lee))))
print(RGB.shape)
plt.figure(figsize=(15, 20))
plt.imshow(RGB)
plt.show()
Sur cette image, vous pouvez voir les taches bleues à Odaiba et Maihama. ** Les points rouges et bleus dans la baie de Tokyo indiquent le navire **, qui montre le navire à la date de prise de vue. Il y a d'autres endroits qui apparaissent en rouge, et il semble que les bâtiments ont disparu ou se sont rassemblés pour devenir une grande chose en raison du développement. On peut dire que la couleur rouge indique également un développement dans ce sens.
Maintenant, enregistrez l'image acquise ici en jpeg avec la date d'observation et les crédits.
import cv2
from PIL import Image, ImageDraw, ImageFont
img_bgr = np.dstack((image_R, np.dstack((image_BG_lee, image_BG_lee))))
new_arr = ((img_bgr - img_bgr.min()) * (1/(img_bgr.max() - img_bgr.min()) * 255)).astype('uint8')
im_rgb = cv2.cvtColor(new_arr, cv2.COLOR_BGR2RGB)
cv2.imwrite(str(object_name) +'.jpg', im_rgb )
img = Image.open(str(object_name) +'.jpg')
plt.figure(figsize=(15, 20))
plt.imshow(img)
plt.show()
En l'état, l'image est blanche dans son ensemble et il est difficile de voir le bleu et le rouge changeants. Maintenant, procédez comme suit pour modifier le contraste. Ceci est votre préféré.
import PIL.Image
from PIL import ImageEnhance
IMAGE_PATH = str(object_name) +'.jpg'
CONTRAST = 2.0
img = PIL.Image.open(IMAGE_PATH)
#Changer le contraste
contrast_converter = ImageEnhance.Contrast(img)
contrast_img = contrast_converter.enhance(CONTRAST)
#Enregistrer l'image
contrast_img.save(str(object_name) +'.jpg')
plt.figure(figsize=(15, 20))
plt.imshow(contrast_img)
plt.show()
Il est plus facile de voir les changements entre le bleu et le rouge qu'auparavant.
Ensuite, la date de prise de vue et le crédit de l'image satellite sont décrits dans cette image. La police utilisée ici est celle qui est ouverte au public gratuitement sur le net. C'est aussi votre préféré.
#Téléchargez et configurez les fichiers de polices
!wget https://osdn.net/dl/mplus-fonts/mplus-TESTFLIGHT-063a.tar.xz
!xz -dc mplus-TESTFLIGHT-*.tar.xz | tar xf -
fontfile = "./mplus-TESTFLIGHT-063a/mplus-1c-bold.ttf"
Maintenant, écrivez les caractères sur l'image.
img = Image.open(str(object_name) +'.jpg')
img = img.convert('RGB')
x = int(img.size[0]/1.21) #Réglage de la position de la description de la date
y = int(img.size[1]/20) #Réglage de la position de la description de la date
fs = int(img.size[0]/70) #Définition de la taille de la police pour la date
obj_draw = ImageDraw.Draw(img)
obj_font = ImageFont.truetype(fontfile, fs)
obj_draw.text((x, y), 'R: '+str(date1), fill=(255, 255, 255), font=obj_font)
obj_draw.text((x, y+1.1*fs), 'G: '+str(date0), fill=(255, 255, 255), font=obj_font)
obj_draw.text((x, y+2.2*fs), 'B: '+str(date0), fill=(255, 255, 255), font=obj_font)
obj_draw.text((img.size[0]/2.0, img.size[1]-y*0.1 - img.size[1]/30 ), 'Contains modified Copernicus Sentinel data (2020)', fill=(255, 255, 255), font=obj_font)
img = img.resize((int(img.size[0] / 2) , int(img.size[1] / 2)))
#img = img.convert('L')
img.save(str(object_name) +'.jpg')
plt.figure(figsize=(15, 20))
plt.imshow(img)
plt.show()
Vous avez maintenant les informations nécessaires sur l'image. Modifiez le code ci-dessus pour ajuster la position, la taille, la couleur, etc. des caractères à votre guise.
Ensuite, créez un ** fichier kmz ** afin que les images acquises ici puissent être visualisées sur un PC local avec une ** application SIG ** (par exemple, ** Google Earth **).
2.6 Créer un fichier kmz à visualiser avec un logiciel SIG
Pour savoir comment créer un fichier kmz, reportez-vous à l'article suivant présenté en détail précédemment.
Je pense qu'il y a plusieurs façons de créer un fichier kmz, mais cette fois j'utiliserai ** simplekml ** fourni par python. Installer et utiliser le module
!pip install simplekml
import simplekml
import zipfile
Ensuite, lisez les informations de latitude / longitude de l'image.
area = pd.DataFrame(AREA,
columns=['longtitude', 'latitude'])
north = area["latitude"].max()
south = area["latitude"].min()
east = area["longtitude"].max()
west = area["longtitude"].min()
Ensuite, créez un fichier kml (informations de latitude / longitude) de la zone d'intérêt ci-dessus.
#Sortie des informations géographiques en kml.
kml = simplekml.Kml()
ground = kml.newgroundoverlay(name=str(object_name))
ground.icon.href = str(object_name) +'.jpg'
ground.latlonbox.north = north
ground.latlonbox.south = south
ground.latlonbox.east = east
ground.latlonbox.west = west
ground.latlonbox.rotation = 0
kml.save(str(object_name)+".kml")
Ensuite, définissez le lien entre le ** fichier kml ** créé et l'image à visualiser, et créez un ** fichier kmz **.
#Génère kmz qui combine des informations géographiques et des images.
with zipfile.ZipFile(str(object_name)+'.kmz', "w", zipfile.ZIP_DEFLATED) as zf:
zf.write(str(object_name)+".kml", arcname=str(object_name)+".kml")
zf.write(str(object_name) +'.jpg', arcname=str(object_name) +'.jpg')
Le fichier kmz est maintenant créé. Vous pouvez consulter le fichier kmz créé par Google Colaboratory. Si vous sélectionnez le côté droit de ce fichier avec la souris, ** download ** sera affiché, alors téléchargez le fichier.
Si ** Google Earth ** est installé sur votre appareil, double-cliquez sur le fichier ** kmz ** téléchargé pour démarrer automatiquement Google Earth et collez l'image créée sur la carte.
Vous pouvez vérifier où l'emplacement est bleu ou rouge sur l'image et désélectionner le fichier kmz pour voir où il se trouve. ** La partie bleue d'Odaiba est-elle l'endroit où mettre le conteneur? Il s'est avéré être **. A part cela, le parking de Tokyo Disneyland est en partie bleu **, ce qui est présumé être le chantier de construction prévu pour la ** nouvelle attraction ** de Disney Sea. (Les photographies aériennes de Google Earth ont un curseur temporel, veuillez donc vérifier la date de prise de vue et les dernières de la partie concernée.)
3. Enfin
Nous avons présenté comment utiliser le moteur Google Earth Engine fourni par Google pour acquérir des images satellite du jour souhaité et voir comment la ville change et où elle change de manière significative. Ici, nous avons évalué principalement la baie de Tokyo et constaté que le développement progressait activement à ** Odaiba et Tokyo Disneyland **. Si vous regardez un autre endroit de la même manière, vous pouvez extraire la ville en développement. De plus, étant donné que le changement d'heure de la situation de la ville peut être observé, les dommages causés par des facteurs de catastrophe autres que le développement peuvent être identifiés par la même méthode.
Personnellement, je pense que la meilleure façon d'utiliser les images satellites est d'extraire ces changements. Les images satellites sont des informations de surface, et il est possible d'obtenir des différences à long terme entre les images passées et présentes et les tendances entre elles (par exemple, le nombre de places de stationnement dans le parking présenté dans l'article précédent). Il y a plusieurs objectifs, mais la méthode est ** l'extraction de changement **. Il est amusant non seulement de parcourir les images satellites, mais aussi de comparer le passé et le présent et de deviner ce qui se passe. J'espère que ce sera l'occasion pour plus de gens de s'intéresser aux images satellites. Si vous avez des commentaires ou des questions, n'hésitez pas à commenter.
Article de référence
J'ai vérifié l'état d'utilisation du parking à partir d'images satellite. J'ai essayé de trouver la tendance du nombre de navires dans la baie de Tokyo à partir d'images satellites. Comment obtenir gratuitement les dernières images satellite.
Kihon de données satellitaires - Compréhension, type, fréquence, résolution, cas d'utilisation ~ Pouvez-vous voir des gens depuis le satellite artificiel? ~ Par satellite, résolution au sol / résumé de l'heure locale ~
Sentinel-1 SAR GRD: C-band Synthetic Aperture Radar Ground Range Detected, log scaling Sentinel-2 MSI: MultiSpectral Instrument, Level-1C
Recommended Posts