lin aerosol

src/CoreRT/CoreRT.jl

@@ -107,6 +107,7 @@ export parameters_from_yaml,                # Getting parameters from a file
        model_from_parameters,               # Converting the parameters to model 
        rt_run,                              # Run the RT code
        default_parameters                   # Set of default parameters
+export lin_model_from_parameters
 
 # Export types to show easily
 export GaussQuadFullSphere, LambertianSurfaceScalar, LambertianSurfaceSpectrum

src/CoreRT/lin_model_from_parameters.jl

@@ -9,7 +9,7 @@ like optical thicknesses, from the input parameters. Produces a vSmartMOM_Model
 default_parameters() = parameters_from_yaml(joinpath(dirname(pathof(vSmartMOM)), "CoreRT", "DefaultParameters.yaml"))
 
 "Take the parameters specified in the vSmartMOM_Parameters struct, and calculate derived attributes into a vSmartMOM_Model" 
-function model_from_parameters(params::vSmartMOM_Parameters)
+function lin_model_from_parameters(params::vSmartMOM_Parameters)
     FT = params.float_type
     #@show FT
     # Number of total bands and aerosols (for convenience)

src/CoreRT/lin_types.jl

@@ -106,45 +106,45 @@ abstract type AbstractLayer end
 "Composite Layer Matrices (`-/+` defined in τ coordinates, i.e. `-`=outgoing, `+`=incoming"
 Base.@kwdef struct linCompositeLayer{FT} <: AbstractLayer 
     "Composite layer Reflectance matrix R (from + -> -)"
-    dR⁻⁺::AbstractArray{FT,3}
+    dR⁻⁺::AbstractArray{FT,4}
     "Composite layer Reflectance matrix R (from - -> +)"
-    dR⁺⁻::AbstractArray{FT,3}
+    dR⁺⁻::AbstractArray{FT,4}
     "Composite layer transmission matrix T (from + -> +)"
-    dT⁺⁺::AbstractArray{FT,3}
+    dT⁺⁺::AbstractArray{FT,4}
     "Composite layer transmission matrix T (from - -> -)"
-    dT⁻⁻::AbstractArray{FT,3}
+    dT⁻⁻::AbstractArray{FT,4}
     "Composite layer source matrix J (in + direction)"
-    dJ₀⁺::AbstractArray{FT,3}
+    dJ₀⁺::AbstractArray{FT,4}
     "Composite layer source matrix J (in - direction)"
-    dJ₀⁻::AbstractArray{FT,3}
+    dJ₀⁻::AbstractArray{FT,4}
 end
 
 "Added (Single) Layer Matrices (`-/+` defined in τ coordinates, i.e. `-`=outgoing, `+`=incoming"
 Base.@kwdef struct linAddedLayer{FT} <: AbstractLayer 
     "Added layer Reflectance matrix R (from + -> -)"
-    dr⁻⁺::AbstractArray{FT,3}
+    dr⁻⁺::AbstractArray{FT,4}
     "Added layer transmission matrix T (from + -> +)"
-    dt⁺⁺::AbstractArray{FT,3}
+    dt⁺⁺::AbstractArray{FT,4}
     "Added layer Reflectance matrix R (from - -> +)"
-    dr⁺⁻::AbstractArray{FT,3}
+    dr⁺⁻::AbstractArray{FT,4}
     "Added layer transmission matrix T (from - -> -)"
-    dt⁻⁻::AbstractArray{FT,3}
+    dt⁻⁻::AbstractArray{FT,4}
     "Added layer source matrix J (in + direction)"
-    dj₀⁺::AbstractArray{FT,3}
+    dj₀⁺::AbstractArray{FT,4}
     "Added layer source matrix J (in - direction)"
-    dj₀⁻::AbstractArray{FT,3}
+    dj₀⁻::AbstractArray{FT,4}
     "Added layer Reflectance matrix R (from + -> -)"
-    dxr⁻⁺::AbstractArray{FT,3}
+    dxr⁻⁺::AbstractArray{FT,4}
     "Added layer transmission matrix T (from + -> +)"
-    dxt⁺⁺::AbstractArray{FT,3}
+    dxt⁺⁺::AbstractArray{FT,4}
     "Added layer Reflectance matrix R (from - -> +)"
-    dxr⁺⁻::AbstractArray{FT,3}
+    dxr⁺⁻::AbstractArray{FT,4}
     "Added layer transmission matrix T (from - -> -)"
-    dxt⁻⁻::AbstractArray{FT,3}
+    dxt⁻⁻::AbstractArray{FT,4}
     "Added layer source matrix J (in + direction)"
-    dxj₀⁺::AbstractArray{FT,3}
+    dxj₀⁺::AbstractArray{FT,4}
     "Added layer source matrix J (in - direction)"
-    dxj₀⁻::AbstractArray{FT,3}
+    dxj₀⁻::AbstractArray{FT,4}
 end
 
 "Composite Layer Matrices (`-/+` defined in τ coordinates, i.e. `-`=outgoing, `+`=incoming"

src/CoreRT/rt_run.jl

@@ -30,7 +30,7 @@ function rt_run(RS_type::AbstractRamanType,
     @unpack quad_points = model
     FT = model.params.float_type
 
-    n_aer = isnothing(model.params.scattering_params) ? 0 : length(model.params.scattering_params.rt_aerosols)
+    Naer = isnothing(model.params.scattering_params) ? 0 : length(model.params.scattering_params.rt_aerosols)
 
     # Also to be changed if more than 1 band is used!!
     # CFRANKEN NEEDS to be changed for concatenated arrays!!
@@ -57,7 +57,7 @@ function rt_run(RS_type::AbstractRamanType,
     dims   = (NquadN,NquadN)              # nxn dims
 
     # Need to check this a bit better in the future!
-    FT_dual = n_aer > 0 ? typeof(model.τ_aer[1][1]) : FT
+    FT_dual = Naer > 0 ? typeof(model.τ_aer[1][1]) : FT
     #FT_dual = FT
     #@show FT_dual
     # Output variables: Reflected and transmitted solar irradiation at TOA and BOA respectively # Might need Dual later!!
@@ -85,12 +85,12 @@ function rt_run(RS_type::AbstractRamanType,
 
     # Create arrays
     @timeit "Creating layers" added_layer         = 
-        make_added_layer(RS_type, FT_dual, arr_type, dims, nSpec)
+        make_added_layer(RS_type, lin, FT_dual, arr_type, dims, nSpec)
     # Just for now, only use noRS here
     @timeit "Creating layers" added_layer_surface = 
-        make_added_layer(RS_type, FT_dual, arr_type, dims, nSpec)
+        make_added_layer(RS_type, lin, FT_dual, arr_type, dims, nSpec)
     @timeit "Creating layers" composite_layer     = 
-        make_composite_layer(RS_type, FT_dual, arr_type, dims, nSpec)
+        make_composite_layer(RS_type, lin, FT_dual, arr_type, dims, nSpec)
     @timeit "Creating arrays" I_static = 
         Diagonal(arr_type(Diagonal{FT}(ones(dims[1]))));
 
@@ -110,7 +110,7 @@ function rt_run(RS_type::AbstractRamanType,
         InelasticScattering.computeRamanZλ!(RS_type, pol_type,Array(qp_μ), m, arr_type)
         # Compute the core layer optical properties:
         @timeit "OpticalProps" layer_opt_props, fScattRayleigh   = 
-            constructCoreOpticalProperties(RS_type,iBand,m,model);
+            constructCoreOpticalProperties(RS_type, iBand, m, model, lin);
         # Determine the scattering interface definitions:
         scattering_interfaces_all, τ_sum_all = 
             extractEffectiveProps(layer_opt_props,quad_points);
@@ -128,7 +128,8 @@ function rt_run(RS_type::AbstractRamanType,
 
             # Expand all layer optical properties to their full dimension:
             @timeit "OpticalProps" layer_opt = 
-                expandOpticalProperties(layer_opt_props[iz], arr_type)
+                expandOpticalProperties(layer_opt_props[iz], 
+                    arr_type)
 
             # Perform Core RT (doubling/elemental/interaction)
             rt_kernel!(RS_type, pol_type, SFI, 
@@ -207,18 +208,23 @@ function rt_run(RS_type::AbstractRamanType,
     #end
 end
 
+function rt_run_test(RS_type::AbstractRamanType, 
+    model::vSmartMOM_Model, lin::vSmartMOM_lin, nparams::Int, 
+    iBand)
+    rt_run(RS_type, model, lin, nparams, iBand)
+end
 
 function rt_run(RS_type::AbstractRamanType, 
-    model::vSmartMOM_Model, lin::vSmartMOM_lin, iBand)
+    model::vSmartMOM_Model, lin::vSmartMOM_lin, nparams::Int, iBand)
     @unpack obs_alt, sza, vza, vaz = model.obs_geom   # Observational geometry properties
     @unpack qp_μ, wt_μ, qp_μN, wt_μN, iμ₀Nstart, μ₀, iμ₀, Nquad = model.quad_points # All quadrature points
     pol_type = model.params.polarization_type
     @unpack max_m = model.params
     @unpack quad_points = model
     FT = model.params.float_type
 
-    n_aer = isnothing(model.params.scattering_params) ? 0 : length(model.params.scattering_params.rt_aerosols)
-
+    Naer = isnothing(model.params.scattering_params) ? 0 : length(model.params.scattering_params.rt_aerosols)
+    Nangles = length(vza)    
     # Also to be changed if more than 1 band is used!!
     # CFRANKEN NEEDS to be changed for concatenated arrays!!
     brdf = model.params.brdf[iBand[1]]
@@ -244,7 +250,7 @@ function rt_run(RS_type::AbstractRamanType,
     dims   = (NquadN,NquadN)              # nxn dims
 
     # Need to check this a bit better in the future!
-    FT_dual = n_aer > 0 ? typeof(model.τ_aer[1][1]) : FT
+    FT_dual = Naer > 0 ? typeof(model.τ_aer[1][1]) : FT
     #FT_dual = FT
     @show FT_dual
     # Output variables: Reflected and transmitted solar irradiation at TOA and BOA respectively # Might need Dual later!!
@@ -309,12 +315,23 @@ function rt_run(RS_type::AbstractRamanType,
 
     # Create arrays
     @timeit "Creating layers" added_layer         = 
-        make_added_layer(RS_type, FT_dual, arr_type, lin, dims, nSpec)
-    # Just for now, only use noRS here
+        make_added_layer(RS_type, FT_dual, arr_type, dims, nSpec)
+    @timeit "Creating layers" lin_added_layer     = 
+        make_lin_added_layer(RS_type, lin, FT, arr_type, dims, nSpec, 
+                    nparams, Int64(3), Int64(4))
+
+        # Just for now, only use noRS here
     @timeit "Creating layers" added_layer_surface = 
-        make_added_layer(RS_type, FT_dual, arr_type, lin, dims, nSpec)
+        make_added_layer(RS_type, FT_dual, arr_type, dims, nSpec)
+    @timeit "Creating layers" lin_added_layer_surface = 
+        make_lin_added_layer(RS_type, lin, FT, arr_type, dims, nSpec, 
+                    nparams, Int64(3), Int64(4))
+
     @timeit "Creating layers" composite_layer     = 
-        make_composite_layer(RS_type, FT_dual, arr_type, lin, dims, nSpec)
+        make_composite_layer(RS_type, FT_dual, arr_type, dims, nSpec)
+    @timeit "Creating layers" lin_composite_layer     = 
+        make_lin_composite_layer(RS_type, lin, FT_dual, arr_type, dims, nSpec, nparams)
+
     @timeit "Creating arrays" I_static = 
         Diagonal(arr_type(Diagonal{FT}(ones(dims[1]))));
 

src/Inelastic/types.jl

@@ -30,6 +30,7 @@ Base.@kwdef mutable struct RRS{FT<:AbstractFloat} <: AbstractRamanType{FT}
     n_Raman::Int
     bandSpecLim = []
     iBand = 1
+    F₀::Array{FT,2} # Solar/Stellar irradiation Stokes vector of size (pol_type.n, nSpec)
 end
 
 """
@@ -57,6 +58,7 @@ Base.@kwdef struct VS_0to1{FT<:AbstractFloat} <: AbstractRamanType{FT}
     dτ₀_λ::FT
     k_Rayl_scatt::FT #σ_Rayl(λ_scatt)/σ_Rayl(λ_incident)
     n_Raman::Int
+    F₀::Array{FT,2} # Solar/Stellar irradiation Stokes vector of size (pol_type.n, nSpec)
 end
 
 """
@@ -83,6 +85,7 @@ Base.@kwdef struct VS_1to0{FT<:AbstractFloat} <: AbstractRamanType{FT}
     dτ₀_λ::FT
     k_Rayl_scatt::FT #σ_Rayl(λ_scatt)/σ_Rayl(λ_incident)
     n_Raman::Int
+    F₀::Array{FT,2} # Solar/Stellar irradiation Stokes vector of size (pol_type.n, nSpec)
     #ramanAtmoProp::RamanAtmosphereProperties
 end
 
@@ -111,13 +114,15 @@ Base.@kwdef mutable struct RVRS{FT<:AbstractFloat} <: AbstractRamanType{FT}
     Z⁺⁺_λ₁λ₀::Array{FT,3} #last dimension -> band index
     τ₀::FT
     n_Raman::Int
+    F₀::Array{FT,2} # Solar/Stellar irradiation Stokes vector of size (pol_type.n, nSpec)
 end
 
 Base.@kwdef mutable struct noRS{FT} <: AbstractRamanType{FT} 
     fscattRayl::Array{FT,1} = [0.0]
     ϖ_Cabannes::Array{FT,1} = [1.0,1.0,1.0] #elastic fraction (Cabannes) of Rayleigh (Cabannes+Raman) scattering
     bandSpecLim = []
     iBand::Array{Int,1} = [1]
+    F₀::Array{FT,2} # Solar/Stellar irradiation Stokes vector of size (pol_type.n, nSpec)
 end
 
 ############################################################
@@ -158,6 +163,7 @@ Base.@kwdef mutable struct RRS_plus{FT<:AbstractFloat} <: AbstractRamanType{FT}
     i_ref::Int
     n_Raman::Int
 
+    F₀::Array{FT,2} # Solar/Stellar irradiation Stokes vector of size (pol_type.n, nSpec)
 end
 
 """
@@ -206,6 +212,7 @@ Base.@kwdef mutable struct VS_0to1_plus{FT<:AbstractFloat} <: AbstractRamanType{
     i_λ₁λ₀_all::Array{Int,1}    = zeros(Int,1)
     i_ref::Int                  = 1
     n_Raman::Int                = 1
+    F₀::Array{FT,2} # Solar/Stellar irradiation Stokes vector of size (pol_type.n, nSpec)
 end
 
 """
@@ -258,13 +265,15 @@ Base.@kwdef mutable struct VS_1to0_plus{FT<:AbstractFloat} <: AbstractRamanType{
     #k_Rayl_scatt::FT #σ_Rayl(λ_scatt)/σ_Rayl(λ_incident)
     n_Raman::Int
     #ramanAtmoProp::RamanAtmosphereProperties
+    F₀::Array{FT,2} # Solar/Stellar irradiation Stokes vector of size (pol_type.n, nSpec)
 end
 
 Base.@kwdef mutable struct noRS_plus{FT} <: AbstractRamanType{FT} 
     fscattRayl::Array{FT,1} = [0.0]
     ϖ_Cabannes::FT = 1.0 #elastic fraction (Cabannes) of Rayleigh (Cabannes+Raman) scattering
     bandSpecLim = []
     iBand::Array{Int,1} = []
+    F₀::Array{FT,2} # Solar/Stellar irradiation Stokes vector of size (pol_type.n, nSpec)
 end
 
 

test/prototyping/Balsamic_OCO2_test.jl

@@ -365,7 +365,7 @@ x = FT[0.2377,     # 1. Legendre Lambertian Albedo, A-band
           690.,    #p₀
           50.]     #σₚ
 
-
+nparams=length(x)
 # Run FW model:
 # @time runner(x);
 #y = oco_sounding.SpectralMeasurement;
@@ -379,7 +379,7 @@ convfct = repeat(h*c*1.e16./λ, Nangles)
 #run_inversion(x,y)
 #run_OE_inversion(x, y, x, Nangles, [length(indices[1]), length(indices[2]), length(indices[3])], convfct);
 Fx = zeros(length(y));
-runner!(Fx,x)
+lin_runner!(Fx,x)
 #Fx = Fx./convfct
 #plot(y, label="Meas") #Ph sec^{-1} m^{-2} sr^{-1} um^{-1}
 #plot!(Fx, label="Mod")# mW/m²-str-cm⁻¹ -> Ph sec^{-1} m^{-2} sr^{-1} um^{-1}