eppaurora.protons.fang2013_protons()   A
last analyzed

Complexity

Conditions 1

Size

Total Lines 49
Code Lines 12

Duplication

Lines 0
Ratio 0 %

Importance

Changes 0
Metric Value
cc 1
eloc 12
nop 5
dl 0
loc 49
rs 9.8
c 0
b 0
f 0
1
# coding: utf-8
2
# Copyright (c) 2020 Stefan Bender
3
#
4
# This file is part of pyeppaurora.
5
# pyeppaurora is free software: you can redistribute it or modify
6
# it under the terms of the GNU General Public License as published
7
# by the Free Software Foundation, version 2.
8
# See accompanying LICENSE file or http://www.gnu.org/licenses/gpl-2.0.html.
9
"""Atmospheric ionization rate parametrizations
10
11
Includes the atmospheric ionization rate parametrization for auroral
12
proton precipitation [1]_.
13
14
.. [1] Fang, X., Lummerzheim, D., and Jackman, C. H. (2013),
15
	Proton impact ionization and a fast calculation method,
16
	J. Geophys. Res. Space Physics, 118, 5369--5378, doi:10.1002/jgra.50484.
17
"""
18
19
import numpy as np
20
from numpy.polynomial.polynomial import polyval
21
22
__all__ = ["fang2013_protons"]
23
24
POLY_F2013 = [
25
	[ 2.55050e+0,  2.69476e-1, -2.58425e-1,  4.43190e-2],
26
	[ 6.39287e-1, -1.85817e-1, -3.15636e-2,  1.01370e-2],
27
	[ 1.63996e+0,  2.43580e-1,  4.29873e-2,  3.77803e-2],
28
	[-2.13479e-1,  1.42464e-1,  1.55840e-2,  1.97407e-3],
29
	[-1.65764e-1,  3.39654e-1, -9.87971e-3,  4.02411e-3],
30
	[-3.59358e-2,  2.50330e-2, -3.29365e-2,  5.08057e-3],
31
	[-6.26528e-1,  1.46865e+0,  2.51853e-1, -4.57132e-2],
32
	[ 1.01384e+0,  5.94301e-2, -3.27839e-2,  3.42688e-3],
33
	[-1.29454e-6, -1.43623e-1,  2.82583e-1,  8.29809e-2],
34
	[-1.18622e-1,  1.79191e-1,  6.49171e-2, -3.99715e-3],
35
	[ 2.94890e+0, -5.75821e-1,  2.48563e-2,  8.31078e-2],
36
	[-1.89515e-1,  3.53452e-2,  7.77964e-2, -4.06034e-3]
37
]
38
39
40
def fang2013_protons(energy, flux, scale_height, rho, pij=None):
41
	"""Proton ionization parametrization by Fang et al., 2013
42
43
	Parametrization for mono-energetic protons as described in [#]_.
44
45
	Parameters
46
	----------
47
	energy: array_like (M,...)
48
		Energy E_0 of the mono-energetic proton beam [keV].
49
	flux: array_like (M,...)
50
		Energy flux Q_0 of the mono-energetic proton beam [keV / cm² / s¹].
51
	scale_height: array_like (N,...)
52
		The atmospheric scale heights [cm].
53
	rho: array_like (N,...)
54
		The atmospheric mass densities [g / cm³], corresponding to the scale heights.
55
	pij: array_like (12, 4), optional
56
		Polynomial coefficents for the proton energy dissipation
57
		per atmospheric depth. Default: `None` (as given in the reference).
58
59
	Returns
60
	-------
61
	en_diss: array_like (M,N)
62
		The dissipated energy profiles [keV].
63
64
	References
65
	----------
66
67
	.. [#] Fang, X., Lummerzheim, D., and Jackman, C. H. (2013),
68
		Proton impact ionization and a fast calculation method,
69
		J. Geophys. Res. Space Physics, 118, 5369--5378, doi:10.1002/jgra.50484.
70
	"""
71
	def _f_y(_c, _y):
72
		# Fang et al., 2008, Eq. (6), Fang et al., 2010 Eq. (4)
73
		# Fang et al., 2013, Eqs. (6), (7)
74
		return (
75
			_c[0] * (_y**_c[1]) * np.exp(-_c[2] * (_y**_c[3])) +
76
			_c[4] * (_y**_c[5]) * np.exp(-_c[6] * (_y**_c[7])) +
77
			_c[8] * (_y**_c[9]) * np.exp(-_c[10] * (_y**_c[11]))
78
		)
79
80
	pij = np.asarray(pij) or np.asarray(POLY_F2013)
81
	# Fang et al., 2013, Eqs. (6), (7)
82
	_cs = np.exp(polyval(np.log(energy), pij.T))
83
	# Fang et al., 2013, Eq. (5)
84
	y = 7.5 / energy * (1e4 * rho * scale_height)**(0.9)
85
	f_y = _f_y(_cs, y)
86
	# Fang et al., 2013, Eq. (3)
87
	en_diss = f_y * flux / scale_height
88
	return en_diss
89