grinsted

Steps

1. use earth engine to export a geotiff of the median modis in higher resolution than you need. https://code.earthengine.google.com/5ce94569601c58a14bcaf28c1852e006
2. download them
3. build a vrt file

gdalbuildvrt g.vrt greenland-modis-median-julyaugust-100m-*.tif

grinsted

import numpy as np
import pandas as pd
import geopandas as gpd
from shapely.geometry import Point, Polygon

zwally = np.loadtxt('grndrainagesystems_ekholm.txt',skiprows=7)

regions = np.unique(zwally[:,0])
for region in regions:
    z = zwally[np.abs(zwally[:,0]-region)<0.001,1:]

grinsted

import rioxarray as rio
from rasterio.enums import Resampling
import xarray as xr
import numpy as np
ds = rio.open_rasterio('BedMachineGreenland-v5.nc')


props_to_save = ["thickness", "surface", "bed", "errbed"]

ds.thickness.rio.write_nodata(0, inplace=True)

grinsted

figure.facecolor: white
figure.dpi: 200
savefig.dpi: 600
figure.figsize: 3.3,2.2
axes.facecolor: white
figure.subplot.bottom: 0.01
figure.subplot.left: 0.01
figure.subplot.right: 0.99
figure.subplot.wspace: 0.2

grinsted

#
# Snippet to empirically determine credible intervals for return period plots
#
# aslak grinsted 2021
#


import numpy as np
from scipy.stats import beta

grinsted

# -*- coding: utf-8 -*-
"""
Created on Thu Feb  6 15:01:40 2020
@author: aslak grinsted
"""


import numpy as np
import matplotlib.pyplot as plt
from numpy import nan

grinsted

import numpy as np
import heapq
import scipy.ndimage


# lebrocq et al.
# https://www.sciencedirect.com/science/article/pii/S0098300406000781?via%3Dihub

def balanceflux(Z, P, n=1):
    """

grinsted

# This file contains code to save 64 samples as a syx file for Novation Circuit
#
# * It reads a single long wave file.
# * splits it into 64 equal chunks
# * crops the sound
# * adds a super fast fadein and "an appropriate fadeout". To remove clicks
# * normalizes every sample
# * Saves the whole thing as a syx file.
#

grinsted

import subprocess
import os
import glob


for filename in glob.iglob(os.path.join(".", "**", "*.nc"), recursive=True):
    if os.stat(filename).st_size < 50:
        continue
    result = subprocess.run(["gdalinfo.exe", filename], stdout=subprocess.PIPE)
    if result.stdout.rfind(b"COMPRESSION=") > 0:

grinsted

from sympy import *
init_printing(use_unicode=False, wrap_line=False, no_global=True)

x = Symbol('x', positive=True)
xmax = Symbol('xmax', positive=True)
sigma = Symbol('sigma', positive=True)
mu = Symbol('mu')

lx = log(x)-log(xmax-x)