eppaurora.electrons.fang2010_maxw_int()   A
last analyzed

Complexity

Conditions 1

Size

Total Lines 46
Code Lines 7

Duplication

Lines 0
Ratio 0 %

Importance

Changes 0
Metric Value
cc 1
eloc 7
nop 7
dl 0
loc 46
rs 10
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 parametrizations for auroral
12
and medium-energy electron precipitation, 100 eV--1 MeV [1]_, [2]_, and [3]_.
13
14
.. [1] Roble and Ridley, Ann. Geophys., 5A(6), 369--382, 1987
15
.. [2] Fang et al., J. Geophys. Res., 113, A09311, 2008
16
.. [3] Fang et al., Geophys. Res. Lett., 37, L22106, 2010
17
"""
18
19
import numpy as np
20
from numpy.polynomial.polynomial import polyval
21
22
from .spectra import pflux_maxwell, ediss_spec_int, ediss_specfun_int
23
24
__all__ = [
25
	"rr1987",
26
	"rr1987_mod",
27
	"fang2008",
28
	"fang2010_mono",
29
	"fang2010_spec_int",
30
	"fang2010_maxw_int",
31
]
32
33
POLY_F2008 = [
34
	[ 3.49979e-1, -6.18200e-2, -4.08124e-2,  1.65414e-2],
35
	[ 5.85425e-1, -5.00793e-2,  5.69309e-2, -4.02491e-3],
36
	[ 1.69692e-1, -2.58981e-2,  1.96822e-2,  1.20505e-3],
37
	[-1.22271e-1, -1.15532e-2,  5.37951e-6,  1.20189e-3],
38
	[ 1.57018,     2.87896e-1, -4.14857e-1,  5.18158e-2],
39
	[ 8.83195e-1,  4.31402e-2, -8.33599e-2,  1.02515e-2],
40
	[ 1.90953,    -4.74704e-2, -1.80200e-1,  2.46652e-2],
41
	[-1.29566,    -2.10952e-1,  2.73106e-1, -2.92752e-2]
42
]
43
44
POLY_F2010 = [
45
	[ 1.24616E+0,  1.45903E+0, -2.42269E-1,  5.95459E-2],
46
	[ 2.23976E+0, -4.22918E-7,  1.36458E-2,  2.53332E-3],
47
	[ 1.41754E+0,  1.44597E-1,  1.70433E-2,  6.39717E-4],
48
	[ 2.48775E-1, -1.50890E-1,  6.30894E-9,  1.23707E-3],
49
	[-4.65119E-1, -1.05081E-1, -8.95701E-2,  1.22450E-2],
50
	[ 3.86019E-1,  1.75430E-3, -7.42960E-4,  4.60881E-4],
51
	[-6.45454E-1,  8.49555E-4, -4.28581E-2, -2.99302E-3],
52
	[ 9.48930E-1,  1.97385E-1, -2.50660E-3, -2.06938E-3]
53
]
54
55
56 View Code Duplication
def rr1987(energy, flux, scale_height, rho):
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
57
	"""Atmospheric electron energy dissipation Roble and Ridley, 1987 [#]_
58
59
	Equations (typo corrected) taken from Fang et al., 2008.
60
61
	Parameters
62
	----------
63
	energy: array_like (M,...)
64
		Characteristic energy E_0 [keV] of the Maxwellian distribution.
65
	flux: array_like (M,...)
66
		Integrated energy flux Q_0 [keV / cm² / s¹]
67
	scale_height: array_like (N,...)
68
		The atmospheric scale heights [cm].
69
	rho: array_like (N,...)
70
		The atmospheric mass density [g / cm³]
71
72
	Returns
73
	-------
74
	en_diss: array_like (M,N)
75
		The dissipated energy profiles [keV].
76
77
	References
78
	----------
79
	.. [#] Roble and Ridley, Ann. Geophys., 5A(6), 369--382, 1987
80
	"""
81
	_c1 = 2.11685
82
	_c2 = 2.97035
83
	_c3 = 2.09710
84
	_c4 = 0.74054
85
	_c5 = 0.58795
86
	_c6 = 1.72746
87
	_c7 = 1.37459
88
	_c8 = 0.93296
89
90
	beta = (rho * scale_height / (4 * 1e-6))**(1 / 1.65)  # RR 1987, p. 371
91
	y = beta / energy  # Corrected in Fang et al. 2008 (4)
92
	f_y = (_c1 * (y**_c2) * np.exp(-_c3 * (y**_c4)) +
93
		_c5 * (y**_c6) * np.exp(-_c7 * (y**_c8)))
94
	# Corrected in Fang et al. 2008 (2)
95
	en_diss = 0.5 * flux / scale_height * f_y
96
	return en_diss
97
98
99 View Code Duplication
def rr1987_mod(energy, flux, scale_height, rho):
0 ignored issues
show
Duplication introduced by
This code seems to be duplicated in your project.
Loading history...
100
	"""Atmospheric electron energy dissipation Roble and Ridley, 1987 [#]_
101
102
	Equations (typo corrected) taken from Fang et al., 2008.
103
	Modified polynomial values to get closer to Fang et al., 2008,
104
	origin unknown.
105
106
	Parameters
107
	----------
108
	energy: array_like (M,...)
109
		Characteristic energy E_0 [keV] of the Maxwellian distribution.
110
	flux: array_like (M,...)
111
		Integrated energy flux Q_0 [keV / cm² / s¹]
112
	scale_height: array_like (N,...)
113
		The atmospheric scale heights [cm].
114
	rho: array_like (N,...)
115
		The atmospheric mass density [g / cm³]
116
117
	Returns
118
	-------
119
	en_diss: array_like (M,N)
120
		The dissipated energy profiles [keV].
121
122
	References
123
	----------
124
	.. [#] Roble and Ridley, Ann. Geophys., 5A(6), 369--382, 1987
125
	"""
126
	# Modified polynomial, origin unknown
127
	_c1 = 3.233
128
	_c2 = 2.56588
129
	_c3 = 2.2541
130
	_c4 = 0.7297198
131
	_c5 = 1.106907
132
	_c6 = 1.71349
133
	_c7 = 1.8835444
134
	_c8 = 0.86472135
135
136
	# Fang et al., 2008, Eq. (4)
137
	y = (rho * scale_height / (4.6 * 1e-6))**(1 / 1.65) / energy
138
	f_y = (_c1 * (y**_c2) * np.exp(-_c3 * (y**_c4)) +
139
		_c5 * (y**_c6) * np.exp(-_c7 * (y**_c8)))
140
	# energy dissipated [keV]
141
	en_diss = 0.5 * flux / scale_height * f_y
142
	return en_diss
143
144
145
def _fang_f_y(_c, _y):
146
	"""Polynomial evaluation helper
147
148
	Fang et al., 2008, Eq. (6), Fang et al., 2010 Eq. (4)
149
	"""
150
	ret = (
151
		_c[0] * (_y**_c[1]) * np.exp(-_c[2] * (_y**_c[3])) +
152
		_c[4] * (_y**_c[5]) * np.exp(-_c[6] * (_y**_c[7]))
153
	)
154
	return ret
155
156
157
def fang2008(energy, flux, scale_height, rho, pij=None):
158
	"""Atmospheric electron energy dissipation from Fang et al., 2008
159
160
	Ionization profile parametrization as derived in Fang et al., 2008 [#]_.
161
162
	Parameters
163
	----------
164
	energy: array_like (M,...)
165
		Characteristic energy E_0 [keV] of the Maxwellian distribution.
166
	flux: array_like (M,...)
167
		Integrated energy flux Q_0 [keV / cm² / s¹]
168
	scale_height: array_like (N,...)
169
		The atmospheric scale height(s) [cm].
170
	rho: array_like (N,...)
171
		The atmospheric densities [g / cm³], corresponding to the scale heights.
172
	pij: array_like (8, 4), optional
173
		Polynomial coefficents for the electron energy dissipation
174
		per atmospheric depth. Default: `None` (as given in the reference).
175
176
	Returns
177
	-------
178
	en_diss: array_like (M,N)
179
		The dissipated energy profiles [keV].
180
181
	References
182
	----------
183
	.. [#] Fang et al., J. Geophys. Res., 113, A09311, 2008, doi: 10.1029/2008JA013384
184
	"""
185
	pij = np.asarray(pij) or np.asarray(POLY_F2008)
186
	# Fang et al., 2008, Eq. (7)
187
	_cs = np.exp(polyval(np.log(energy), pij.T))
188
	# Fang et al., 2008, Eq. (4)
189
	y = (rho * scale_height / (4e-6))**(1 / 1.65) / energy
190
	f_y = _fang_f_y(_cs, y)
191
	# Fang et al., 2008, Eq. (2)
192
	en_diss = 0.5 * f_y * flux / scale_height
193
	return en_diss
194
195
196
def fang2010_mono(energy, flux, scale_height, rho, pij=None):
197
	r"""Atmospheric electron energy dissipation from Fang et al., 2010
198
199
	Parametrization for mono-energetic electrons [#]_.
200
201
	Parameters
202
	----------
203
	energy: array_like (M,...)
204
		Energy E_0 of the mono-energetic electron beam [keV].
205
	flux: array_like (M,...)
206
		Energy flux Q_0 of the mono-energetic electron beam [keV / cm² / s¹].
207
	scale_height: array_like (N,...)
208
		The atmospheric scale heights [cm].
209
	rho: array_like (N,...)
210
		The atmospheric mass densities [g / cm³], corresponding to the scale heights.
211
	pij: array_like (8, 4), optional
212
		Polynomial coefficents for the electron energy dissipation
213
		per atmospheric depth. Default: `None` (as given in the reference).
214
215
	Returns
216
	-------
217
	en_diss: array_like (M,N)
218
		The dissipated energy profiles [keV].
219
220
	References
221
	----------
222
	.. [#] Fang et al., Geophys. Res. Lett., 37, L22106, 2010, doi: 10.1029/2010GL045406
223
	"""
224
	pij = np.asarray(pij) or np.asarray(POLY_F2010)
225
	# Fang et al., 2010, Eq. (5)
226
	_cs = np.exp(polyval(np.log(energy), pij.T))
227
	# Fang et al., 2010, Eq. (1)
228
	y = 2. / energy * (rho * scale_height / (6e-6))**(0.7)
229
	f_y = _fang_f_y(_cs, y)
230
	# Fang et al., 2008, Eq. (2)
231
	en_diss = f_y * flux / scale_height
232
	return en_diss
233
234
235
def fang2010_spec_int(ens, dfluxes, scale_height, rho, pij=None, axis=-1):
236
	r"""Integrate over a given energy spectrum
237
238
	Integrates over the mono-energetic parametrization `q` from [#]_
239
	using the given differential particle spectrum `phi`:
240
241
	:math:`\int_\text{spec} \phi(E) q(E, Q) E \text{d}E`
242
243
	Parameters
244
	----------
245
	ens: array_like (M,...)
246
		Central (bin) energies of the spectrum
247
	dfluxes: array_like (M,...)
248
		Differential particle fluxes in the given bins
249
	scale_height: array_like (N,...)
250
		The atmospheric scale heights
251
	rho: array_like (N,...)
252
		The atmospheric densities, corresponding to the
253
		scale heights.
254
	pij: array_like (8, 4), optional
255
		Polynomial coefficents for the electron energy dissipation
256
		per atmospheric depth. Default: `None` (as given in the reference).
257
	axis: int, optional
258
		The axis to use for integration, default: -1 (last axis).
259
260
	Returns
261
	-------
262
	en_diss: array_like (N)
263
		The dissipated energy profiles [keV].
264
265
	References
266
	----------
267
	.. [#] Fang et al., Geophys. Res. Lett., 37, L22106, 2010, doi: 10.1029/2010GL045406
268
269
	See Also
270
	--------
271
	fang2010_mono, ediss_spec_int
272
	"""
273
	return ediss_spec_int(
274
		ens, dfluxes, scale_height, rho, fang2010_mono,
275
		axis=axis,
276
		func_kws=dict(pij=pij),
277
	)
278
279
280
def fang2010_maxw_int(energy, flux, scale_height, rho, bounds=(0.1, 300.), nstep=128, pij=None):
281
	"""Integrate Fang et al., 2010 over a Maxwellian spectrum
282
283
	Integrates the mono-energetic parametrization from Fang et al., 2010 [#]_
284
	over a Maxwellian spectrum with characteristic energy `energy` and
285
	total energy flux `flux`.
286
287
	Parameters
288
	----------
289
	energy: float or array_like (M,...)
290
		Characteristic energy E_0 [keV] of the Maxwellian distribution.
291
	flux: float or array_like (M,...)
292
		Integrated energy flux Q_0 [keV / cm² / s¹]
293
	scale_height: float or array_like (N,...)
294
		The atmospheric scale heights [cm].
295
	rho: float or array_like (N,...)
296
		The atmospheric mass density [g / cm³]
297
	bounds: tuple, optional
298
		(min, max) [keV] of the integration range to integrate the Maxwellian.
299
		Make sure that this is appropriate to encompass the spectrum.
300
		Default: (0.1, 300.)
301
	nsteps: int, optional
302
		Number of integration steps, default: 128.
303
	pij: array_like (8, 4), optional
304
		Polynomial coefficents for the electron energy dissipation
305
		per atmospheric depth. Default: `None` (as given in the reference).
306
307
	Returns
308
	-------
309
	en_diss: array_like (M,N)
310
		The dissipated energy profiles [keV].
311
312
	References
313
	----------
314
	.. [#] Fang et al., Geophys. Res. Lett., 37, L22106, 2010, doi: 10.1029/2010GL045406
315
316
	See Also
317
	--------
318
	fang2010_mono, fang2010_specfun_int, pflux_maxwell
319
	"""
320
	return ediss_specfun_int(
321
		energy, flux, scale_height, rho, fang2010_mono,
322
		ediss_kws=dict(pij=pij),
323
		bounds=bounds,
324
		nstep=nstep,
325
		spec_fun=pflux_maxwell,
326
	)
327