Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Question about POROSITY and DENST #39

Closed
pikaqiu2002 opened this issue Dec 18, 2024 · 4 comments
Closed

Question about POROSITY and DENST #39

pikaqiu2002 opened this issue Dec 18, 2024 · 4 comments

Comments

@pikaqiu2002
Copy link

def bulk_ice_density(temperature, salinity, porosity):
    """
    Computes bulk density of sea ice (in kg m :sup:`-3`), when considering the influence from  brine, solid salts, and
    air bubbles in the ice. Formulation from Cox & Weeks (1983): Equations for determining the gas and brine volumes in sea ice samples,
    J Glac. Developed for temperatures between -2--30°C. For higher temperatures (>2°C) is used the formulation from
    Lepparanta & Manninen (1988): The brine and gas content of sea ice with attention to low salinities and high temperatures.

    :param temperature: Temperature in K
    :param salinity: salinity in kg/kg (see PSU constant in smrt module)
    :param porosity: Fractional volume of air inclusions (0..1)
    :returns: Density of ice mixture in kg m :sup:`-3`
    """

    Tc = temperature - FREEZING_POINT

    if Tc > -2.0:
        alpha = np.array([-4.1221e-2, -18.407, 5.8402e-1, 2.1454e-1])
        beta = np.array([9.0312e-2, -1.6111e-2, 1.2291e-4, 1.3603e-4])

    elif Tc >= -22.9:
        alpha = np.array([-4.732, -22.45, -6.397e-1, -1.074e-2])
        beta = np.array([8.903e-2, -1.763e-2, -5.33e-4, -8.801e-6])

    else:
        alpha = np.array([9.899e3, 1.309e3, 55.27, 7.160e-1])
        beta = np.array([8.547, 1.089, 4.518e-2, 5.819e-4])

    F1 = np.polyval(alpha[::-1], Tc)
    F2 = np.polyval(beta[::-1], Tc)

    # Density of pure ice, C&W p. 311:
    rho_ice = 0.917 - 1.403e-4 * Tc  # in g/cm^3

    # Density of mixture:
    rho = (1. - porosity) * (rho_ice * F1 / (F1 - rho_ice * salinity * PSU**-1 * F2)) * 1e3  # in kg/m3 (eq. 15, C&W, 1983)

    if rho < 0:
        raise SMRTError("Ice density may not be negative.")

    return rho

1996cf29601a8e8ac615ff3cc0d9478
Good afternoon all!

I've set the density parameters and they are all positive, but I still get the error. I wonder if when the code is run, it is still running with the value of porosity=0 instead of using the density value I entered.

Thank you very much!

@JulienMeloche
Copy link
Contributor

SMRT will calculate the porosity using the function bulk_ice_density with porosity = 0 since you entered the density and no porosity was set (see line 505 in inputs.make_medium.py). I'm pretty sure its either the temperature or salinity you set. Did you input the temperature in Kelvin or Celsius? Check also the unit of your salinity and PSU.

Its looks like whatever is coming out of the density of mixture equation is negative... Which is based on the temperature and salinity your entered.

@pikaqiu2002
Copy link
Author

SMRT will calculate the porosity using the function bulk_ice_density with porosity = 0 since you entered the density and no porosity was set (see line 505 in inputs.make_medium.py). I'm pretty sure its either the temperature or salinity you set. Did you input the temperature in Kelvin or Celsius? Check also the unit of your salinity and PSU.

Its looks like whatever is coming out of the density of mixture equation is negative... Which is based on the temperature and salinity your entered.

Thank you for your answer! I will check temperature and salinity.
I have another question. Can SMRT output the penetration depth of electromagnetic waves in the sea ice in a specific active sensor mode?
Thank you very much!

@JulienMeloche
Copy link
Contributor

Yes, using the result instance, you can calculate the extinction coefficient (ke = ks+ka) from the scattering (ks) and absorption coefficient (ka). The penetration depth (or the e-folding depth) is defined by 1/ke.

for exemple:

result = model.run(sensor, snowpack)

you can calculate ke and the penetration depth with :

ke = result.other_data['ks'] + result.other_data['ka']
d_p = ke**-1

This calculation is only dependent on the frequency of the sensor, not the type (active, passive or altimeter).

You could also use the optical depth. More info in the doc:
https://smrt.readthedocs.io/en/master/smrt.core.html#smrt.core.result.Result.optical_depth

@pikaqiu2002
Copy link
Author

Yes, using the result instance, you can calculate the extinction coefficient (ke = ks+ka) from the scattering (ks) and absorption coefficient (ka). The penetration depth (or the e-folding depth) is defined by 1/ke.

for exemple:

result = model.run(sensor, snowpack)

you can calculate ke and the penetration depth with :

ke = result.other_data['ks'] + result.other_data['ka']
d_p = ke**-1

This calculation is only dependent on the frequency of the sensor, not the type (active, passive or altimeter).

You could also use the optical depth. More info in the doc: https://smrt.readthedocs.io/en/master/smrt.core.html#smrt.core.result.Result.optical_depth

Thank you for your answer!

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants