Check Mie

RemoteSensingTools · Oct 29, 2024 · aae70e0 · aae70e0
1 parent 87467a1
commit aae70e0
Show file tree

Hide file tree

Showing 35 changed files with 6,712 additions and 161 deletions.
diff --git a/Project.toml b/Project.toml
@@ -4,13 +4,18 @@ authors = ["Rupesh Jeyaram <[email protected]> and contributors"]
 version = "0.5.0"
 
 [deps]
+AWS = "fbe9abb3-538b-5e4e-ba9e-bc94f4f92ebc"
+AWSS3 = "1c724243-ef5b-51ab-93f4-b0a88ac62a95"
+ArgParse = "c7e460c6-2fb9-53a9-8c5b-16f535851c63"
+Base64 = "2a0f44e3-6c83-55bd-87e4-b1978d98bd5f"
 BenchmarkTools = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf"
 CUDA = "052768ef-5323-5732-b1bb-66c8b64840ba"
 CUDAKernels = "72cfdca4-0801-4ab0-bf6a-d52aa10adc57"
 CairoMakie = "13f3f980-e62b-5c42-98c6-ff1f3baf88f0"
 ClimaParams = "5c42b081-d73a-476f-9059-fd94b934656c"
 ColorSchemes = "35d6a980-a343-548e-a6ea-1d62b119f2f4"
 Colors = "5ae59095-9a9b-59fe-a467-6f913c188581"
+DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"
 DataInterpolations = "82cc6244-b520-54b8-b5a6-8a565e85f1d0"
 Decimals = "abce61dc-4473-55a0-ba07-351d65e31d42"
 DelimitedFiles = "8bb1440f-4735-579b-a4ab-409b98df4dab"
@@ -21,15 +26,20 @@ DocStringExtensions = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae"
 FastGaussQuadrature = "442a2c76-b920-505d-bb47-c5924d526838"
 Formatting = "59287772-0a20-5a39-b81b-1366585eb4c0"
 ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210"
+GLMakie = "e9467ef8-e4e7-5192-8a1a-b1aee30e663a"
 GR = "28b8d3ca-fb5f-59d9-8090-bfdbd6d07a71"
+GeoMakie = "db073c08-6b98-4ee5-b6a4-5efafb3259c6"
+GeometryBasics = "5c1252a2-5f33-56bf-86c9-59e7332b4326"
 GridLayoutBase = "3955a311-db13-416c-9275-1d80ed98e5e9"
+HDF5 = "f67ccb44-e63f-5c2f-98bd-6dc0ccc4ba2f"
 HTTP = "cd3eb016-35fb-5094-929b-558a96fad6f3"
 ImageFiltering = "6a3955dd-da59-5b1f-98d4-e7296123deb5"
 Insolation = "e98cc03f-d57e-4e3c-b70c-8d51efe9e0d8"
 InstrumentOperator = "9e589c1b-9e01-4e00-831a-aa39ce86e3ef"
 Interpolations = "a98d9a8b-a2ab-59e6-89dd-64a1c18fca59"
 JLD2 = "033835bb-8acc-5ee8-8aae-3f567f8a3819"
 JSON = "682c06a0-de6a-54ab-a142-c8b1cf79cde6"
+JSON3 = "0f8b85d8-7281-11e9-16c2-39a750bddbf1"
 KernelAbstractions = "63c18a36-062a-441e-b654-da1e3ab1ce7c"
 LaTeXStrings = "b964fa9f-0449-5b57-a5c2-d3ea65f4040f"
 LegendrePolynomials = "3db4a2ba-fc88-11e8-3e01-49c72059a882"
@@ -39,15 +49,19 @@ Missings = "e1d29d7a-bbdc-5cf2-9ac0-f12de2c33e28"
 NCDatasets = "85f8d34a-cbdd-5861-8df4-14fed0d494ab"
 NNlib = "872c559c-99b0-510c-b3b7-b6c96a88d5cd"
 NaNMath = "77ba4419-2d1f-58cd-9bb1-8ffee604a2e3"
+NearestNeighbors = "b8a86587-4115-5ab1-83bc-aa920d37bbce"
 NetCDF = "30363a11-5582-574a-97bb-aa9a979735b9"
 NetCDF_jll = "7243133f-43d8-5620-bbf4-c2c921802cf3"
 OffsetArrays = "6fe1bfb0-de20-5000-8ca7-80f57d26f881"
 Parameters = "d96e819e-fc66-5662-9728-84c9c7592b0a"
 Pkg = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f"
+PlotUtils = "995b91a9-d308-5afd-9ec6-746e21dbc043"
+PlotlyBase = "a03496cd-edff-5a9b-9e67-9cda94a718b5"
 PlotlyJS = "f0f68f2c-4968-5e81-91da-67840de0976a"
 Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
 Polynomials = "f27b6e38-b328-58d1-80ce-0feddd5e7a45"
 ProgressMeter = "92933f4c-e287-5a05-a399-4b506db050ca"
+PyCall = "438e738f-606a-5dbb-bf0a-cddfbfd45ab0"
 PyPlot = "d330b81b-6aea-500a-939a-2ce795aea3ee"
 Revise = "295af30f-e4ad-537b-8983-00126c2a3abe"
 SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b"

diff --git a/SSA.png b/SSA.png
diff --git a/extXS.png b/extXS.png
diff --git a/imag_refind.png b/imag_refind.png
diff --git a/local_alb.png b/local_alb.png
diff --git a/nonvegbiome_SIF.png b/nonvegbiome_SIF.png
diff --git a/real_refind.png b/real_refind.png
diff --git a/rho1.png b/rho1.png
diff --git a/scattXS.png b/scattXS.png
diff --git a/scatter_plot.png b/scatter_plot.png
diff --git a/shettle_1.dat~ b/shettle_1.dat~
diff --git a/shettle_2.dat~ b/shettle_2.dat~
diff --git a/shettle_3.dat~ b/shettle_3.dat~
diff --git a/shettle_4.dat~ b/shettle_4.dat~
diff --git a/shettle_5.dat~ b/shettle_5.dat~
diff --git a/src/CoreRT/atmo_prof.jl b/src/CoreRT/atmo_prof.jl
@@ -26,7 +26,7 @@ function compute_atmos_profile_fields(T, p_half::AbstractArray, q, vmr; g₀=9.8
     # Now actually compute the layer VCDs
     for i = 1:n_layers 
         Δp = p_half[i + 1] - p_half[i]
-        vmr_h2o[i] = dry_mass/(dry_mass-wet_mass*(1-1/q[i]))
+        vmr_h2o[i] = (dry_mass/wet_mass)*q[i]/(1-q[i])
         vmr_dry = 1 - vmr_h2o[i]
         M  = vmr_dry * dry_mass + vmr_h2o[i] * wet_mass
         vcd = Nₐ * Δp / (M  * g₀ * 100^2) * 100
@@ -444,4 +444,4 @@ function compute_absorption_profile!(τ_abs::Array{FT,2},
         τ_abs[:,iz] += Array(absorption_cross_section(absorption_model, grid, p, T)) * profile.vcd_dry[iz] * vmr_curr
     end
 
-end
+end
diff --git a/src/CoreRT/model_from_parameters.jl b/src/CoreRT/model_from_parameters.jl
@@ -138,7 +138,7 @@ function model_from_parameters(params::vSmartMOM_Parameters)
     greek_cabannes = Vector{vSmartMOM.Scattering.GreekCoefs{Float64}}()
     greek_rayleigh = Vector{vSmartMOM.Scattering.GreekCoefs{Float64}}()
     ϖ_Cabannes = zeros(n_bands)
-    # Remove rayleight for testing:
+    # Remove rayleigh for testing:
     τ_rayl = [zeros(params.float_type,length(params.spec_bands[i]), length(profile.T)) for i=1:n_bands];
     #τ_rayl = [zeros(params.float_type,1,length(profile.T)) for i=1:n_bands];
 
@@ -186,33 +186,61 @@ function model_from_parameters(params::vSmartMOM_Parameters)
         =#
 
                                 #@show τ_rayl[i_band]
+
         # If no absorption, continue to next band
-        isnothing(params.absorption_params) && continue
-        #@show "hello"
-        #@show params.absorption_params.molecules[i_band]
-        #@show length(params.absorption_params.molecules[i_band])
-        # Loop over all molecules in this band, obtain profile for each, and add them up
-        for molec_i in 1:length(params.absorption_params.molecules[i_band])
-
-            # This can be precomputed as well later in my mind, providing an absorption_model or an interpolation_model!
-            if isempty(params.absorption_params.luts)
-                # Obtain hitran data for this molecule
-                @timeit "Read HITRAN" hitran_data = read_hitran(artifact(params.absorption_params.molecules[i_band][molec_i]), iso=1)
+        (isnothing(params.absorption_params) && isnothing(params.q)) && continue
 
-                println("Computing profile for $(params.absorption_params.molecules[i_band][molec_i]) with vmr $(profile.vmr[params.absorption_params.molecules[i_band][molec_i]]) for band #$(i_band)")
+        if !isnothing(params.q)
+            # Obtain hitran data for this H₂O
+            @timeit "Read HITRAN" hitran_data = read_hitran(artifact("H2O"), iso=1)
+
+            println("Computing profile for water vapor in band #$(i_band)")
                 # Create absorption model with parameters beforehand now:
                 absorption_model = make_hitran_model(hitran_data, 
-                    params.absorption_params.broadening_function, 
-                    wing_cutoff = params.absorption_params.wing_cutoff, 
-                    CEF = params.absorption_params.CEF, 
+                    vSmartMOM.Absorption.Voigt(), 
+                    wing_cutoff = 150, 
+                    CEF = vSmartMOM.Absorption.HumlicekWeidemann32SDErrorFunction(), 
                     architecture = params.architecture, 
                     vmr = 0);#mean(profile.vmr[params.absorption_params.molecules[i_band][molec_i]]))
                 # Calculate absorption profile
 
-                @timeit "Absorption Coeff"  compute_absorption_profile!(τ_abs[i_band], absorption_model, params.spec_bands[i_band],profile.vmr[params.absorption_params.molecules[i_band][molec_i]], profile);
-            # Use LUT directly
-            else
-                compute_absorption_profile!(τ_abs[i_band], params.absorption_params.luts[i_band][molec_i], params.spec_bands[i_band],profile.vmr[params.absorption_params.molecules[i_band][molec_i]], profile);
+                @timeit "Absorption Coeff"  compute_absorption_profile!(τ_abs[i_band], 
+                    absorption_model, 
+                    params.spec_bands[i_band],
+                    profile.vmr_h2o, 
+                    profile);
+        end
+        #@show "hello"
+        #@show params.absorption_params.molecules[i_band]
+        #@show length(params.absorption_params.molecules[i_band])
+        # Loop over all molecules in this band, obtain profile for each, and add them up
+        if !isnothing(params.absorption_params)
+            for molec_i in 1:length(params.absorption_params.molecules[i_band])
+
+                # This can be precomputed as well later in my mind, providing an absorption_model or an interpolation_model!
+                if isempty(params.absorption_params.luts)
+                    # Obtain hitran data for this molecule
+                    @timeit "Read HITRAN" hitran_data = read_hitran(artifact(params.absorption_params.molecules[i_band][molec_i]), iso=1)
+
+                    println("Computing profile for $(params.absorption_params.molecules[i_band][molec_i]) with vmr $(profile.vmr[params.absorption_params.molecules[i_band][molec_i]]) for band #$(i_band)")
+                    # Create absorption model with parameters beforehand now:
+                    absorption_model = make_hitran_model(hitran_data, 
+                        params.absorption_params.broadening_function, 
+                        wing_cutoff = params.absorption_params.wing_cutoff, 
+                        CEF = params.absorption_params.CEF, 
+                        architecture = params.architecture, 
+                        vmr = 0);#mean(profile.vmr[params.absorption_params.molecules[i_band][molec_i]]))
+                    # Calculate absorption profile
+
+                    @timeit "Absorption Coeff"  compute_absorption_profile!(τ_abs[i_band], 
+                        absorption_model, 
+                            params.spec_bands[i_band],
+                                profile.vmr[params.absorption_params.molecules[i_band][molec_i]], 
+                                profile);
+                # Use LUT directly
+                else
+                    compute_absorption_profile!(τ_abs[i_band], params.absorption_params.luts[i_band][molec_i], params.spec_bands[i_band],profile.vmr[params.absorption_params.molecules[i_band][molec_i]], profile);
+                end
             end
         end
     end

diff --git a/src/CoreRT/show_utils.jl b/src/CoreRT/show_utils.jl
@@ -78,8 +78,8 @@ function Base.show(io::IO,::MIME"text/plain", x::vSmartMOM_Parameters)
             #println(io, "\t\t  σ: $(curr_rt_aerosol.aerosol.size_distribution.σ) μm; geometric standard deviation")
             println(io, "\t\t  nᵣ: $(curr_rt_aerosol.aerosol.nᵣ)")
             println(io, "\t\t  nᵢ: $(curr_rt_aerosol.aerosol.nᵢ)")
-            println(io, "\t\t  p₀: $(curr_rt_aerosol.p₀) Pa") 
-            println(io, "\t\t  σp: $(curr_rt_aerosol.σp) Pa")
+            println(io, "\t\t  z₀: $(curr_rt_aerosol.z₀) km") 
+            println(io, "\t\t  σ₀: $(curr_rt_aerosol.σ₀) ")
         end
         println(io, "\tr_max: $(x.scattering_params.r_max) μm")
         println(io, "\tnquad_radius: $(x.scattering_params.nquad_radius)")

diff --git a/src/SIF_emission/ficus_refl.dat~ b/src/SIF_emission/ficus_refl.dat~