Commit c3e73768 authored by Tilmann Sager's avatar Tilmann Sager
Browse files

Added generation of plots

parent e5c0d91a
from glob import glob
import geopandas as gpd
import matplotlib.pyplot as plt
import pandas as pd
import seaborn as sns
from shapely.geometry import Point
from constants import Columns as col
sns.set_theme(style="whitegrid")
"""
PLOTS
"""
def boxplot(data, filename):
plot = sns.boxplot(x=data)
fig = plot.get_figure()
fig.savefig(filename)
def histogram(img, filename):
plot = sns.histplot(img.ravel(), bins=256)
fig = plot.get_figure()
fig.savefig(filename)
pass
def density_map(title, flight_gdf, contrail_gdf):
world_map = gpd.read_file('../contraildetect/data/maps/ne_10m_admin_0_countries/ne_10m_admin_0_countries.shp')
fig, ax = plt.subplots(figsize=(16, 10))
world_map.to_crs(epsg=4326).plot(ax=ax, color='lightgrey')
flight_gdf.plot(ax=ax, cmap='Reds',
legend=True, legend_kwds={'shrink': 0.3},
markersize=1)
contrail_gdf.plot(ax=ax, cmap='Greens', legend=False, markersize=1)
ax.set_title(title)
plt.savefig(title + '.png')
def get_flight_gdf():
all_flights = pd.DataFrame()
for flight_csv in glob('../contraildetect/data/flights/processed/*.csv'):
all_flights = all_flights.append(pd.read_csv(flight_csv))
crs = {'init': 'EPSG:4326'}
geometry = [Point(xy) for xy in zip(all_flights['LONGITUDE'], all_flights['LATITUDE'])]
return gpd.GeoDataFrame(crs=crs, geometry=geometry)
def get_contrail_gdf(df):
crs = {'init': 'EPSG:4326'}
geometry = [Point(xy) for xy in df[col.line_coord].values.tolist()]
return gpd.GeoDataFrame(crs=crs, geometry=geometry)
def _normalize_column(df: pd.DataFrame, column_name: str):
return (df[column_name] - df[column_name].min()) / (
df[column_name].max() - df[column_name].min())
def test_boxplot():
data_csv = '/home/noam/src/python/contraildetect/data/output/2021-06-08/results.csv'
data = pd.read_csv(data_csv)
flight_counts = data[data['number_flights'] > 0]['number_flights']
boxplot(flight_counts, 'flight_counts.png')
def _flight_pixel_relation(row: pd.Series):
result = 0
if row[col.px_count_norm] != 0 and row[col.flight_count_norm]:
result = row[col.px_count_norm] / row[col.flight_count_norm]
return result
def test_contrail_map(df):
flights = get_flight_gdf()
contrails = get_contrail_gdf(df)
density_map('contrail_flights', flights, contrails)
def flight_pixel_relation(df: pd.DataFrame):
df[col.flight_count_norm] = _normalize_column(df, col.flight_count)
df[col.px_count_norm] = _normalize_column(df, col.px_count)
df[col.rel_flight_px] = df.apply(_flight_pixel_relation, axis=1)
return df
......@@ -13,16 +13,19 @@ class Columns:
dim = 'dimension'
flight_raw = 'flight_raw'
flight_processed = 'flight_processed'
flight_count = 'number_of_flights'
flight_count = 'flight_count'
flight_count_norm = 'flight_count_norm'
bands = 'bands'
cda = 'cda'
cda_t = 'cda_t'
lines = 'lines'
line_arr = 'line_arr'
line_count = 'line_count'
line_coord = 'line_center_coordinate'
line_coord = 'line_coordinates'
line_center = 'line_center_coordinates'
px_dens = 'pixel_density'
px_count = 'pixel_count'
px_count_norm = 'pixel_count_norm'
rel_flight_px = 'flight_to_pixel_relation'
......@@ -34,7 +37,11 @@ columns_to_save = [Columns.hdf,
Columns.line_count,
Columns.px_dens,
Columns.px_count,
Columns.rel_flight_px]
Columns.rel_flight_px,
Columns.west,
Columns.east,
Columns.north,
Columns.south]
DTYPES = {
'FLIGHT_ID': np.int,
......
import os
from glob import glob
import pandas as pd
from tqdm import tqdm
import analysis
import config
import plots
from constants import columns_to_save
from stages import prepare, transformation, detection, flight_filtering, segmentation, stats
from util import check_and_repair_csv
......@@ -21,6 +23,7 @@ input_files = glob(config.product_dir() + '/*.hdf')
# TODO: multiprocessing
# TODO: implement opt-out/in stages
# TODO: improve logging
# TODO: filter time per day
print('Creating job')
job_df = prepare.create_job(input_files)
......@@ -67,10 +70,28 @@ ASSEMBLING INFORMATION
print('Assembling information')
job_df = job_df.progress_apply(stats.collect, axis=1)
"""
ANALYSIS
"""
print('Running analysis')
job_df = analysis.flight_pixel_relation(job_df)
"""
RESULTS
"""
print('Saving results')
job_df.to_csv(os.path.join(config.output_dir(), 'output.csv'), columns=columns_to_save, index=False)
analysis.test_contrail_map(job_df)
"""
PLOTS
"""
print('Creating plots')
all_flights = pd.DataFrame()
for flight_csv in glob(config.processed_flights() + '/*.csv'):
all_flights = all_flights.append(pd.read_csv(flight_csv))
# plots.boxplot(job_df)
# plots.flights(all_flights)
# plots.granules(job_df)
plots.contrails(job_df)
plots.density_map(job_df, all_flights)
from stages import transformation, detection, segmentation, assembling
from stages import transformation, detection, segmentation
def run(granule):
......
import geopandas as gpd
import geoplot as gplt
import matplotlib.pyplot as plt
import pandas as pd
import seaborn as sns
from pyproj import CRS
from shapely.geometry import Point, Polygon, LineString
from constants import Columns as col
sns.set_theme(style="whitegrid")
world_map = gpd.read_file('../contraildetect/data/maps/ne_10m_admin_0_countries/ne_10m_admin_0_countries.shp')
crs = CRS.from_epsg(4326)
default_include = ['flights', 'granules', 'contrails']
alpha = 0.1
def boxplot(df: pd.DataFrame, filename: str = 'flight_counts.png'):
flight_counts = df[df[col.flight_count] > 0][col.flight_count]
plot = sns.boxplot(x=flight_counts)
fig = plot.get_figure()
fig.savefig(filename)
def histogram(img, filename):
plot = sns.histplot(img.ravel(), bins=256)
fig = plot.get_figure()
fig.savefig(filename)
pass
def _flight_gdf(df: pd.DataFrame, visualize_as: str):
geometry = []
if visualize_as == 'pt':
geometry = [Point(xy) for xy in zip(df['LONGITUDE'], df['LATITUDE'])]
if visualize_as == 'line':
geometry = []
return gpd.GeoDataFrame(crs=crs, geometry=geometry)
def flights(flight_df: pd.DataFrame, visualize_as: str = 'pt', filename: str = 'flight_density.png', cmap='Reds'):
fig, ax = plt.subplots(figsize=(16, 10))
world_map.to_crs(epsg=4326).plot(ax=ax, color='lightgrey')
gdf = _flight_gdf(flight_df, visualize_as)
gplt.kdeplot(gdf, ax=ax, cmap=cmap, shade=True, legend=False, alpha=0.5, markersize=1)
# gdf.plot(ax=ax, cmap=cmap, legend=False, alpha=0.5, markersize=1)
plt.savefig(filename)
def _contrail_gdf(df: pd.DataFrame, visualize_as: str):
geometry = []
for index, row in df.iterrows():
for line in row[col.line_coord]:
if visualize_as == 'line':
geometry.append(LineString([Point(line[0]), Point(line[1])]))
if visualize_as == 'pt':
geometry.append(Point(line[0]))
geometry.append(Point(line[1]))
return gpd.GeoDataFrame(crs=crs, geometry=geometry)
def contrails(df: pd.DataFrame, visualize_as: str = 'line', filename: str = 'contrails.png', cmap: str = 'Blues'):
fig, ax = plt.subplots(figsize=(16, 10))
world_map.to_crs(epsg=4326).plot(ax=ax, color='lightgrey')
gdf = _contrail_gdf(df, visualize_as)
gdf.plot(ax=ax, cmap=cmap, legend=False, alpha=alpha)
plt.savefig(filename)
def _granule_gdf(df: pd.DataFrame, visualize_as: str):
points_lst = []
for index, row in df.iterrows():
lat_min = min(row[col.north], row[col.south])
lat_max = max(row[col.north], row[col.south])
lon_min = min(row[col.west], row[col.east])
lon_max = max(row[col.west], row[col.east])
points_lst.append(
[Point(lon_min, lat_min), Point(lon_min, lat_max), Point(lon_max, lat_max), Point(lon_max, lat_min)])
geometry = [Polygon(points) for points in points_lst]
return gpd.GeoDataFrame(df, crs=crs, geometry=geometry)
def granules(df: pd.DataFrame, visualize_as: str = 'poly', filename: str = 'granules.png', cmap: str = 'Blues'):
fig, ax = plt.subplots(figsize=(16, 10))
world_map.to_crs(epsg=4326).plot(ax=ax, color='lightgrey')
gdf = _granule_gdf(df, visualize_as)
gdf.plot(ax=ax, cmap=cmap, legend=False, alpha=alpha)
plt.savefig(filename)
def density_map(df: pd.DataFrame, flight_df: pd.DataFrame,
include: [str] = default_include, title: str = 'density_map'):
fig, ax = plt.subplots(figsize=(16, 10))
world_map.to_crs(epsg=4326).plot(ax=ax, color='lightgrey')
if 'granules' in include:
granule_gdf = _granule_gdf(df, 'TODO')
granule_gdf.plot(column=col.rel_flight_px, ax=ax, cmap="Blues", legend=False, alpha=alpha)
if 'flights' in include:
flight_gdf = _flight_gdf(flight_df, 'pt')
gplt.kdeplot(flight_gdf, ax=ax, cmap='Reds', shade=True, legend=False, alpha=0.5, markersize=1)
if 'contrails' in include:
contrail_gdf = _contrail_gdf(df, 'line')
contrail_gdf.plot(ax=ax, cmap='Greens', legend=False, alpha=alpha)
ax.set_title(title)
plt.savefig(title + '.png')
"""
north-south: latitude
west-east: longitude
"""
import numpy as np
import pandas as pd
from constants import Columns as col
from stages.mapping import map_point_to_coordinates
......@@ -31,21 +32,27 @@ def relation_flight_number_pixel_density(number_flights, pixel_density):
return 0
def line_coordinates(lines, dim, lat, lon):
def line_center(lines, dim, lat, lon):
lines_coordinates = []
if lines:
for line in lines:
line_x = (line[0][0], line[1][0])
line_y = (line[0][1], line[1][1])
line_center = (max(line_x) - min(line_x), max(line_y) - min(line_y))
lines_coordinates.append(map_point_to_coordinates(line_center, dim, lon, lat))
line_cent = (max(line_x) - min(line_x), max(line_y) - min(line_y))
lines_coordinates.append(map_point_to_coordinates(line_cent, dim, lon, lat))
return lines_coordinates
def map_lines_to_coordinates():
return 0
def map_lines_to_coordinates(lines, dim, lat, lon):
lines_coordinates = []
if lines:
for line in lines:
line_start = map_point_to_coordinates(line[0], dim, lon, lat)
line_end = map_point_to_coordinates(line[1], dim, lon, lat)
lines_coordinates.append((line_start, line_end))
return lines_coordinates
# TODO
......@@ -55,8 +62,11 @@ def collect(granule):
granule[col.px_dens] = relative_pixel_density(granule)
granule[col.rel_flight_px] = relation_flight_number_pixel_density(granule[col.flight_count],
granule[col.px_dens])
granule[col.line_coord] = line_coordinates(granule[col.lines], granule[col.dim],
(granule[col.north], granule[col.south]),
(granule[col.west], granule[col.east]))
granule[col.line_center] = line_center(granule[col.lines], granule[col.dim],
(granule[col.north], granule[col.south]),
(granule[col.west], granule[col.east]))
granule[col.line_coord] = map_lines_to_coordinates(granule[col.lines], granule[col.dim],
(granule[col.north], granule[col.south]),
(granule[col.west], granule[col.east]))
return granule
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment