smrt-model / smrt

Snow Microwave Radiative Transfert model to compute thermal emission and backscatter from snowpack
Other
51 stars 20 forks source link

Using the model for non snow packs #3

Closed Anthony-CamembearGames closed 2 years ago

Anthony-CamembearGames commented 5 years ago

Hello,

Not sure if this is the best place to ask.

Is it possible to use the model for layers composed of other things than snow ?

My current approach to this is to use make_generic_stack as follows:

# For a case with no scattering

from smrt import make_model, core, sensor_list
import numpy as np
from smrt.core.globalconstants import C_SPEED
from smrt.inputs.make_medium import make_generic_stack

sensor = sensor_list.passive(562e9, 40)
l = 2
nl = l//2 # // Forces integer division
thickness = np.array([0.1, 0.1]*nl)
thickness[-1] = 100  # last one is semi-infinit
temperature = np.array([150.0, 150.0]*nl)
permitivity = np.array([3+0.1*1j, 2.5+0.08*1j])

k0 = 2 * np.pi * sensor.frequency / C_SPEED
ka = 2 * k0 * np.sqrt(permitivity).imag
ks = 0

genericpack = make_generic_stack(thickness=thickness, 
                          temperature=temperature, 
                          ks=ks, 
                          ka=ka, 
                          effective_permittivity=permitivity)

m = make_model("prescribed_kskaeps", "dort")
result = m.run(sensor, genericpack)

print(np.mean(result.Tb().values))

Now if I want to add scattering to the model what would be the best way to proceed ? I assume it is not simply a case of changing the ka and ks values as such (for Rayleigh scattering):

lmda = C_SPEED / sensor.frequency

radius = 1e-4
eps = permitivity 
e0 = 1 # vacuum

f = 0.6 # the faction of scatteres such that 1-f = porosity

k0 = 2 * np.pi / lmda
ks = f * 2 * abs((eps - e0) / (eps + 2 * e0))**2 * radius**3 * k0**4
ka = f * 9 * k0 * eps.imag / e0 * abs(e0 / (eps + 2 * e0))**2

Is this a correct approach ?

Thank you very much for your help,

ghislainp commented 5 years ago

Dear Anthony,

your way to do it is ok, it should work, but another (more 'SMRT') way is as follow:

make_snowpack(thickness, "independent_sphere", ice_permittivity_model=eps_of_your_scatterer_material)

where eps_of_your_scatterer_materiel is the permittivity of your particular material (ice is used by default but the keyword ice_permittivity_model allows overwriting the default)

Then

make_model("rayleigh", "dort")

This is how SMRT is supposed to be used for rayleigh scattering. Of course after trying that this works, a more elegant approach is to create a small wrapper function "make_mymedium" that takes meaningful input arguments for your domain, and call behind the scene make_snowpack or even directly Layer() constructor. The recommended way to improve SMRT is indeed to add as many functions as necessary similarly to make_medium.py calling Layer() to build the layers. I even recommend to not add this function in the smrt structure but instead outside in your own module so you can benefit from future SMRT improvments (and soon some dramatic speed increase will come!). At least the idea is to never-never touch files in the core/ directory. This should not be necessary. But of course this is only recommendation if you want to keep long-term compatibility with SMRT.

Best regards,

Ghislain

On 2/14/19 8:31 PM, Anthony Lethuillier wrote:

Hello,

Not sure if this is the best place to ask.

Is it possible to use the model for layers composed of other things than snow ?

My current approach to this is to use |make_generic_stack| as follows:

For a case with no scattering

|from smrt import make_model, core, sensor_list import numpy as np from smrt.core.globalconstants import C_SPEED from smrt.inputs.make_medium import make_generic_stack sensor = sensor_list.passive(562e9, 40) l = 2 nl = l//2 # // Forces integer division thickness = np.array([0.1, 0.1]nl) thickness[-1] = 100 # last one is semi-infinit temperature = np.array([150.0, 150.0]nl) permitivity = np.array([3+0.11j, 2.5+0.081j]) k0 = 2 np.pi sensor.frequency / C_SPEED ka = 2 k0 np.sqrt(permitivity).imag ks = 0 genericpack = make_generic_stack(thickness=thickness, temperature=temperature, ks=ks, ka=ka, effective_permittivity=permitivity) m = make_model("prescribed_kskaeps", "dort") result = m.run(sensor, genericpack) print(np.mean(result.Tb().values)) |

Now if I want to add scattering to the model what would be the best way to proceed ? I assume it is not simply a case of changing the |ka| and |ks| values as such (for Rayleigh scattering):

|lmda = C_SPEED / sensor.frequency radius = 1e-4 eps = permitivity e0 = 1 # vacuum f = 0.6 # the faction of scatteres such that 1-f = porosity k0 = 2 np.pi / lmda ks = f 2 * abs((eps - e0) / (eps + 2

  • e0))2 * radius*3 k04 ka = f 9 k0 eps.imag / e0 abs(e0 / (eps + 2 * e0))**2 |

Thank you very much for your help,

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/smrt-model/smrt/issues/3, or mute the thread https://github.com/notifications/unsubscribe-auth/AKEwHhuqEANyKpj8_vr5fzRH8-QunpLCks5vNbmbgaJpZM4a8SXi.

Anthony-CamembearGames commented 5 years ago

Thanks for your answer,

I was testing the method you suggested and I have an additional question,

If I use the make_snowpack function it requires the density and radius parameters to be fixed. If I use this snowpack with a non-scattering model the results do not change with the radius but they do change with the density. Additionally the results are different than when I use the make_generic_stack function.

See the following code for an example of this:

from smrt import make_snowpack, make_model, core, sensor_list
import numpy as np
from smrt.core.globalconstants import C_SPEED
from smrt.inputs.make_medium import make_generic_stack

sensor = sensor_list.passive(562e9, 40)
l = 2
nl = l//2 # // Forces integer division
thickness = np.array([0.1, 0.1]*nl)
thickness[-1] = 100  # last one is semi-infinit
temperature = np.array([250.0, 250.0]*nl)
permitivity = np.array([3+0.1*1j, 3+0.1*1j])
density = np.array([200.0, 200.0]*nl)
radius = np.array([1e-5, 1e-5]*nl)

def eps_of_your_scatterer_material(frequency, temp):
    return 3+0.1*1j

genericpack = make_snowpack(thickness=thickness, 
              temperature=temperature, 
              density = density,
              radius = radius,
              microstructure_model="independent_sphere", 
              ice_permittivity_model=eps_of_your_scatterer_material)

m = make_model("nonscattering", "dort")
result = m.run(sensor, genericpack)

print('Snowpack Model:')
print(np.mean(result.Tb().values))

k0 = 2 * np.pi * sensor.frequency / C_SPEED
ka = 2 * k0 * np.sqrt(permitivity).imag
ks = 0

genericpack = make_generic_stack(thickness=thickness, 
                          temperature=temperature, 
                          ks=ks, 
                          ka=ka, 
                          effective_permittivity=permitivity)

m = make_model("prescribed_kskaeps", "dort")
result = m.run(sensor, genericpack)

print('Generic stack Model:')
print(np.mean(result.Tb().values))

I am correct in thinking that the density is only used to calculate the fraction of vacuum and solid ? So in my case I would have to define the DENSITY_OF_ICE variable to be equal to the bulk density of my material ?

Regards,

Anthony

Anthony-CamembearGames commented 5 years ago

Ok I find that if I change density = np.array([DENSITY_OF_ICE, DENSITY_OF_ICE]*nl) the two temperatures become equal.