irfpy.util.fluxtools

Tools for flux (differential, velocity distribution)

Code author: Yoshifumi Futaana

This module provide tools for flux manipulations (differential flux, count rate, velocity distribution function)

Containers of differential flux and count rate in a spherical coordinates are implemented in irfpy.util.energyspectrum. The velocity distribution functions are implemented in irfpy.util.vdf module.

irfpy.util.energyspectrum.DifferentialFluxGrid(...)

Representation of the differential flux on grids in polar coordinates.

irfpy.util.energyspectrum.CountRateGrid(cps, ...)

Represent the count rate (counts per sec) on grids in polar frame.

irfpy.util.vdf.VelocityDistributionGrid(vdf, ...)

Represents the velocity distribution function on grids in polar frame

vdf2dflux(vdfun)

Conversion from velocity distribution function to differential flux.

dflux2vdf(dflx)

Conversion from differential flux to velocity distribution function.

cnt2dflx_simple_convfact(count_rate, gfactor)

Conversion from count rate to differential flux using g-factor.

Code author: Yoshifumi Futaana

irfpy.util.fluxtools.vdf2dflux(vdfun)[source]

Conversion from velocity distribution function to differential flux.

Parameters:

vdfun (irfpy.util.vdf.VelocityDistributionGrid) – Velocity distribution function on polar grid.

Returns:

Differential flux object.

Return type:

irfpy.util.energyspectrum.DifferentialFluxGrid

In MKSA, the conversion from VDF to Dflx is

\[J = \frac{2 f E}{m^2}\]

Limitation: So far, the energy/theta/phi_grid_size is limited to default (see irfpy.util.energyspectrum.DifferentialFluxGrid. Only the central value is used for conversion.

>>> vgrid = np.linspace(0.01, 1000e3, 96)   # 96 velocity space
>>> tgrid = np.linspace(30, 150, 32)     # 32 polar steps
>>> pgrid = np.linspace(-180, 180, 65)[:-1]  # 64 azimuth steps
>>> maxwell = vdf.maxwellVelocityDistributionGrid(1e6, [-400e3, 0, 0], 80e3, vgrid, np.deg2rad(tgrid), np.deg2rad(pgrid))
>>> maxwell_dflx = vdf2dflux(maxwell)
>>> print('%.4f' % maxwell_dflx.density())
0.9999
>>> v = maxwell_dflx.velocity_vector()
>>> print('%.2f %.2f %.2f' % (v[0], v[1], v[2]))
-400.02 -0.00 0.00
irfpy.util.fluxtools.dflux2vdf(dflx)[source]

Conversion from differential flux to velocity distribution function.

Parameters:

dflx – Differential flux object in irfpy.util.energyspectrum.DifferentialFluxGrid object.

Returns:

Velocity distribution function in irfpy.util.vdf.VelocityDistributionGrid class.

>>> vgrid = np.linspace(0.01, 1000e3, 96)   # 96 velocity space
>>> tgrid = np.linspace(30, 150, 32)     # 32 polar steps
>>> pgrid = np.linspace(-180, 180, 65)[:-1]  # 64 azimuth steps
>>> maxwell = vdf.maxwellVelocityDistributionGrid(1e6, [-400e3, 0, 0], 80e3, vgrid, np.deg2rad(tgrid), np.deg2rad(pgrid))
>>> maxwell_dflx = vdf2dflux(maxwell)
>>> maxwell_vdf = dflux2vdf(maxwell_dflx)   # Back again.
>>> print('%.5f' % ((maxwell.data - maxwell_vdf.data) ** 2).max())
0.00000
irfpy.util.fluxtools.cnt2dflx_simple_convfact(count_rate, gfactor)[source]

Conversion from count rate to differential flux using g-factor.

Conversion from count rate to differential flux using g-factor (table) is implemented.

\[J = \frac{C}{G \times E}\]

where C is the count rate [/s], G is the geometric factor including efficiency [cm2 sr eV/eV], and E is the energy [eV]. The obtained differential flux is J, with the unit of [/cm2 sr eV s].

The usage is identical to dflx2cnt_simple_convfact(), which provides inverse function.

Parameters:

  • count_rate (irfpy.util.energyspectrum.CountRateGrid) – Calculated count rate.

  • gfactor – Geometric factor, in the same shape as count_rate.data (iene, ipol, iaz), or broadcast-able. If gfactor depends on the energy, and gfactor has the shape of (iene), which is not broadcastable. In such case, you may need to make the gfactor array by using irfpy.util.nptools.broadcast_r(). See the sample below.

Returns:

Differential flux object.

Return type:

irfpy.util.energyspectrum.DifferentialFluxGrid

irfpy.util.fluxtools.dflx2cnt_simple_convfact(diff_flux, gfactor)[source]

Differential flux is converted to count rate.

With a simple conversion of “by definition”.

\[C = J \times G \times E\]

where J is the differential flux [/cm2 s sr eV], G is the geometric factor including efficiency [cm2 sr eV/eV], E is the energy [eV] (not energy range). The returned count rate will be [/s].

Parameters:

  • diff_flux (irfpy.util.energyspectrum.DifferentialFluxGrid) – Differential flux object.

  • gfactor – Geometric factor, in the same shape as dflx.data (iene, ipol, iaz), or broadcast-able. If gfactor depends on the energy, and gfactor has the shape of (iene), which is not broadcastable. In such case, you may need to make the gfactor array by using irfpy.util.nptools.broadcast_r(). See the sample below.

Returns:

Calculated count rate.

Return type:

irfpy.util.energyspectrum.CountRateGrid

For comprehensive example, I first prepare the energy array for a “virtual” sensor:

>>> elist = [10, 20, 40, 80, 160, 320, 640, 1280, 2560, 5120]    # 10 energy steps

Then, angle information of the sensor are defined.

>>> tlist = np.radians([15, 45, 75, 105, 135, 165])    # Polar angle, 6 channel
>>> plist = np.radians([0, 45, 90, 135, 180, 225, 270, 315])   # Azimuth, 8 channels

As an exercise, I prepare the Maxwell distribution

>>> n = 1e6   # 1 /cc
>>> vb = [200e3, 200e3, 0]  # 200 km/s vx and vy
>>> vth = 70e3  # 70 km/s thermal velocity
>>> diff_flux = energyspectrum.maxwellDifferentialFluxGrid(n, vb, vth, elist, tlist, plist)
>>> print('{:.2f}'.format( diff_flux.data.max()))
68765.73

Now the diff_flux is an object of irfpy.util.energyspectrum.DifferentialFluxGrid

Giving the g-factor of 1e-5 cm2 sr eV/eV, you may get the count rate object (irfpy.util.energyspectrum.CountRateGrid object).

>>> cps = dflx2cnt_simple_convfact(diff_flux, 1e-5)
>>> print('{:.1f}'.format(cps.cps.max()))
440.1

If the g-factor has energy dependence, and you have g-factor array, you may need to make a g-factor array to be (nEne, nPol, nAz) shape.

>>> gfactor = np.array([1, 2, 3, 4, 5, 5, 4, 3, 2, 1]) * 1e-5
>>> gfactor_3d = nptools.broadcast_r(gfactor, [10, 6, 8])
>>> print(gfactor_3d.shape)
(10, 6, 8)
>>> cps = dflx2cnt_simple_convfact(diff_flux, gfactor_3d)
>>> print('{:.1f}'.format(cps.cps.max()))
1760.4

The inversion function is also prepared as cnt2dflx_simple_convfact().

>>> dflx2 = cnt2dflx_simple_convfact(cps, gfactor_3d)

The results should be (almost) identical.

>>> diff = dflx2.dflux - diff_flux.dflux
>>> print('{:.10f}'.format((diff ** 2).max()))
0.0000000000
irfpy.util.fluxtools.fit_differential_flux_by_maxwell(diff_flux, log=False, initial_params=None, method='L-BFGS-B')[source]

Fit the energy spectrum by Maxwellian using least-square method.

Parameters:
Author:

Yoshifumi Futaana

>>> estep = [10, 20, 40, 80, 160, 320, 640, 1280, 2560, 5120, 10240, 20480]
>>> tstep = np.deg2rad([22.5, 45, 67.5, 90, 112.5, 135, 157.5])
>>> pstep = np.deg2rad([0, 30, 60, 90, 120, 150, 180, 210, 240, 270, 300, 330])
>>> maxwell_flux = energyspectrum.maxwellDifferentialFluxGrid(
...                     10e6, [400e3, 100e3, -50e3], 200e3,
...                     estep, tstep, pstep)
>>> fit = fit_differential_flux_by_maxwell(maxwell_flux)
>>> print('{:.1f} {:.1f} {:.1f} {:.1f} {:.1f}'.format(fit[0], fit[1], fit[2], fit[3], fit[4]))
10.0 400.0 100.0 -50.0 200.0