Source code for hydrac.model.scattering.hp

"""
High-pass models (:mod:`hydrac.model.scattering.hp`)
====================================================

.. autoclass:: Hp
   :members:
   :private-members:

"""

import numpy as np


[docs]class Hp(object): """High-pass models class. This class deals with high-pass model formulations of exact models in the Rayleigh and geometric regions. They are easy to compute due to their parametric description. Parameters ---------- obj : object Scattering object from :class: `hydrac.model.scattering.scattering.Scattering` f : float (numpy.array) Instrument frequencies a : float (numpy.array) Size vector """ def __init__(self, obj, f, a): self.obj = obj self.obj.scattering.atmp = a self.obj.scattering.theta = self.obj.particle.param.ori self.obj.scattering.katmp = 2 * np.pi * np.outer(f, a) / self.obj.cw self.f = f if hasattr(self.obj.particle, 'irregular_tag') is False: self.obj.particle.irregular_tag = input('Canonical (CAN) or ' + 'irregular (IRR) shape ?') g = self.obj.particle.param.rhos / self.obj.rho h = self.obj.particle.param.c / self.obj.cw self.alpha_pis = (1 - g * h**2) / (3 * g * h**2) + (1 - g) / (1+2 * g) self.alpha_pic = (1 - g * h**2) / (2 * g * h**2) + (1 - g) / (1 + g) self.R = (g * h - 1) / (g * h + 1) self.obj.scattering.f_inftmp,\ self.obj.scattering.c_inftmp = self.hp_func()
[docs] def hp_func(self): """ This method computes the Stanton (1989) models for all types of particles. This model is only valid for angles close to backscattering, ie for an angle between the incident pressure field and the principal axis of the particle :math:`\\theta \\approx -\\pi/2` (this is only true for cylinders or prolate spheroids). Returns ------- f_inf : float (numpy.array) Model form function :math:`f_{\infty}` c_inf : float (numpy.array) Model total scattering cross-secton :math:`\\chi_{\infty}` """ a = self.obj.scattering.atmp ka = self.obj.scattering.katmp try: L = a * self.obj.particle.param.L_a kL = ka * self.obj.particle.param.L_a s = np.sin(kL * np.cos(self.obj.scattering.theta * np.pi / 180)) /\ (kL * np.cos(self.obj.scattering.theta * np.pi / 180)) except AttributeError: pass Ka = ka * np.sin(self.obj.scattering.theta * np.pi / 180) if self.obj.particle.irregular_tag == 'CAN': G = 1 F = 1 if self.obj.particle.shape == 'SPH': if self.obj.particle.irregular_tag == 'IRR': if self.obj.particle.nature == 'MIN': F = 15 * ka**(-1.9) G = 1 elif self.obj.particle.nature == 'ORG': F = 40 * ka**(-0.4) G = 1 - 0.8 * np.exp(-2.5 * (ka - 2.25)**2) if self.obj.particle.nature == 'GAS': F = 1 + 0.5 * ka**(-0.75) G = 1 + 85 * np.exp(-500000 * (ka - 0.0135)**2) sig = a**2 * ka**4 * self.alpha_pis**2 * G / (1+4 * ka**4 * self.alpha_pis**2 / (self.R**2 * F)) elif self.obj.particle.shape == 'PRO': if self.obj.particle.irregular_tag == 'IRR': if self.obj.particle.nature == 'MIN': F = 1.8 * ka**(-0.4) G = 1 elif self.obj.particle.nature == 'ORG': F = 2.5 * ka**1.65 G = 1 - 0.8 * np.exp(-2.5 * (ka - 2.3)**2) if self.obj.particle.nature == 'GAS': F = 1 G = 1 + 6 * np.exp(-700000 * (ka - 0.0045)**2) sig = (1 / 9) * L**2 * ka**4 * self.alpha_pic**2 * G /\ (1+(16/9) * ka**4 * self.alpha_pic**2 / (self.R**2 * F)) elif self.obj.particle.shape == 'CYL': if self.obj.particle.irregular_tag == 'IRR': if self.obj.particle.nature == 'MIN': F = 3.5 * (ka**(-1)) G = 1 elif self.obj.particle.nature == 'ORG': F = 3.0 * ka**0.65 G = 1 - 0.8 * np.exp(-2.5 * (ka - 2.0)**2) if self.obj.particle.nature == 'GAS': F = 1 + 0.35 * ka**(-0.8) G = 1 + 20 * np.exp(-900000 * (ka-0.0045)**2) sig = (1/4) * L**2 * Ka**4 * self.alpha_pic**2 * s**2 * G /\ (1 + np.pi * Ka**3 * self.alpha_pic**2 / (self.R**2 * F)) elif self.obj.particle.shape == 'bCYL': H = 0.5 + 0.5 * self.obj.particle.param.rho_L *\ np.sin(self.obj.particle.param.rho_L**(-1)) if self.obj.particle.irregular_tag == 'IRR': if self.obj.particle.nature == 'MIN': F = 2.5 * (ka**(-1)) G = 1 elif self.obj.particle.nature == 'ORG': F = 3.0 * ka**0.65 G = 1 - 0.8 * np.exp(-2.5 * (ka - 2.0)**2) if self.obj.particle.nature == 'GAS': F = 1 + 0.1 * ka**(-0.75) G = 1 + 7 * np.exp(-800000 * (ka - 0.0045)**2) sig = (1/4) * L**2 * ka**4 * self.alpha_pic**2 * H**2 * G /\ (1 + L**2 * ka**4 * self.alpha_pic**2 * H**2 / (self.obj.particle.param.rho_L * L * a * self.R**2 * F)) if self.obj.particle.shape in ['CYL', 'bCYL']: return 4 * np.sqrt(sig) / np.tile(np.atleast_2d(L), (sig.shape[0], 1)),\ np.zeros(np.shape(sig)) else: return 2 * np.sqrt(sig) / np.tile(np.atleast_2d(a), (sig.shape[0], 1)),\ np.zeros(np.shape(sig))