import numpy as np
import math
import matplotlib.pyplot as plt
from scipy.stats import norm
from scipy.stats import crystalball

# Good function to model "lossy" tails on low side
# Crystall ball: f(x,beta,m)  (beta>0 and m>1)
# for x>-beta it is a gaussian exp(-x**2/)
# for x<-beta it is a power law (B-x)**m
# (so beta is effectively in units of sigma)

# Compare Crystall Ball with Gaussian
mean  = 100
sigma = 5
beta  = 1
m     = 2

x  = np.linspace(max(0,mean-6*sigma), mean+4*sigma, 1000) # plot range
ax = plt.subplot(111)
ax.plot(x, crystalball.pdf(x, beta, m, scale=sigma, loc=mean), label="Crystall Ball")
ax.plot(x, norm.pdf(x, scale=sigma, loc=mean),  label="Gaussian")
ax.grid()
ax.legend()
_ = ax.set_ylim(bottom=0)