diff --git a/examples/7_SLM_Fraunhofer.jl b/examples/7_SLM_Fraunhofer.jl new file mode 100644 index 0000000..80d09bb --- /dev/null +++ b/examples/7_SLM_Fraunhofer.jl @@ -0,0 +1,215 @@ +### A Pluto.jl notebook ### +# v0.19.41 + +using Markdown +using InteractiveUtils + +# This Pluto notebook uses @bind for interactivity. When running this notebook outside of Pluto, the following 'mock version' of @bind gives bound variables a default value (instead of an error). +macro bind(def, element) + quote + local iv = try Base.loaded_modules[Base.PkgId(Base.UUID("6e696c72-6542-2067-7265-42206c756150"), "AbstractPlutoDingetjes")].Bonds.initial_value catch; b -> missing; end + local el = $(esc(element)) + global $(esc(def)) = Core.applicable(Base.get, el) ? Base.get(el) : iv(el) + el + end +end + +# ╔═╡ d551f366-0213-11ef-176b-7fcf5b47315e +begin + using Pkg + Pkg.activate(".") + Pkg.instantiate() + using Revise +end + +# ╔═╡ 9e6d909f-06de-46e7-929a-61bfd67bea26 +using Zygote, WaveOpticsPropagation, FFTW, ImageShow, FiniteDifferences, CUDA, Optim, IndexFunArrays, FourierTools, TestImages + +# ╔═╡ ef79854f-0ae9-41aa-8f92-aae8d7023d71 +using Plots + +# ╔═╡ f9a13d72-f904-467c-8346-a9ef272693ff +using PlutoUI + +# ╔═╡ 20c5faa9-5694-42fd-b1bb-0c5093f6c93c +begin + # use CUDA if functional + use_CUDA = Ref(true && CUDA.functional()) + var"@mytime" = use_CUDA[] ? CUDA.var"@time" : var"@time" + + togoc(x) = use_CUDA[] ? CuArray(x) : x +end + +# ╔═╡ b8cf117e-d084-4d42-8022-03fd6ac4c417 +TableOfContents() + +# ╔═╡ 94ea2891-db37-477e-b535-1a6d402d4af2 + + +# ╔═╡ a971769f-1215-4422-a209-fad710be16dd +sz = (1024, 1024) + +# ╔═╡ bdf3b801-28e0-41cf-9dc4-92fe8fdaca0e +beam = togoc(0im .+ gaussian(Float32, sz, sigma=100)) .* togoc(cispi.(conv((normal(Float32, (sz), sigma=5)), 3 .* rand(Float32, (sz))))); + +# ╔═╡ 36b798c9-feee-43b7-b9ec-4f0558a44d2c +simshow(Array(beam)) + +# ╔═╡ ad4ef631-841d-412e-a09a-4db60fa7e4a9 +z = 100f-3 + +# ╔═╡ 15011619-f26a-4858-b0a0-d5dbcad49824 +λ = 633f-9 + +# ╔═╡ 48048a82-5bbf-4f16-a7d9-15993c83f905 +L = 10f-3 + +# ╔═╡ 3f791571-2a36-4479-a853-c6ed4e8a270d +md"# Forward model" + +# ╔═╡ e52cfcd8-eec8-4a8b-b953-f1c6ef2be959 +N = 3 + +# ╔═╡ 42c00ce2-d1e4-4f84-84be-a7e4cd6551ce +begin + diffuser = cispi.(conv((normal(Float32, (sz), sigma=1)), 10 .* rand(Float32, (sz..., N)))) + diffuser_c = togoc(diffuser)# .* box(sz, (120, 150)) + + diffuser_c[:, :, 1] .= cis.(0.7 .* rr2(sz, scale=0.05)) .+ 0.1f0 .* diffuser_c[:, :, 1] + diffuser_c[:, :, 2] .= cis.(0.7 .*rr2(sz, scale=0.05)) .+ 0.1f0 .* diffuser_c[:, :, 2] + diffuser_c[:, :, 3] .= cis.(0.7 .*rr2(sz, scale=0.05)) .+ 0.1f0 .* diffuser_c[:, :, 3] + + #diffuser_c[:, :, 3] .= 1 + +# diffuser_c[:, :, 4] .*= cis.(0.3 .* xx((sz))) +# diffuser_c[:, :, 5] .*= cis.(.- 0.3 .* xx((sz))) +# diffuser_c[:, :, 6] .*= cis.(0.3 .* yy((sz))) +# diffuser_c[:, :, 7] .*= cis.(.- 0.3 .* yy((sz))) +end + +# ╔═╡ d058c9ba-82ec-4f1e-8938-eb9e1e95b713 +fraunhofer = Fraunhofer(beam .* diffuser_c, z, λ, L) + +# ╔═╡ 5d21e0bc-ea6e-4b2c-8977-d1761d1b3c65 +fraunhofer(beam .* diffuser_c) ./ 0.00195312 ≈ fftshift(fft(ifftshift(beam .* diffuser_c, (1,2)), (1,2)), (1,2)) + +# ╔═╡ df72f725-8464-46b0-a291-8503ea4e563f +size(diffuser_c) + +# ╔═╡ 90d78713-0d4d-48b2-b590-a79d4e7c12d5 +simshow(Array(diffuser_c[:, :, 1])) + +# ╔═╡ f9c4d5ca-a35a-462a-b5f4-4120012f5a1b +fwd(x) = begin + return abs2.(fraunhofer(x .* diffuser_c)) +end + +# ╔═╡ d84457d3-30ae-4711-b47f-7d2a7376a237 +md"# Optimizer" + +# ╔═╡ d219be8b-2d22-4fbe-b022-217edfeec5dc +function make_fg!(fwd, measured, loss=:L2) + L2 = let measured=measured, fwd=fwd + function L2(x) + return sum(abs2, fwd(x) .- measured) + end + end + + f = let + if loss == :L2 + L2 + end + end + + g! = let f=f + function g!(G, rec) + if !isnothing(G) + return G .= Zygote.gradient(f, rec)[1] + end + end + end + return f, g! +end + +# ╔═╡ 692e0263-8273-45a6-b3e6-c30bbeaeb649 +measurement_c = fwd(beam); + +# ╔═╡ 5c550775-a614-44dc-a81d-638448d23067 +@bind iz Slider(1:N) + +# ╔═╡ 56a03baf-3873-43d5-b8c1-8569de6026b0 +simshow(Array(measurement_c)[:, :, iz], γ=1) + +# ╔═╡ fc9bb298-5b46-4585-9ea5-4df59405d18b +fraunhofer.params.Lp + +# ╔═╡ de16936c-7774-4aab-9702-f0cf4fa660fd +md"# Optimize" + +# ╔═╡ ee5260f3-d3a4-4bf2-8b54-31723c3c516f +rec0 = togoc(ones(ComplexF32, sz)); + +# ╔═╡ d90bff7e-150d-407d-a510-385fd5a32df3 +f, g! = make_fg!(fwd, measurement_c) + +# ╔═╡ eadee3cb-1302-4f1a-a21c-17265d9a2e6b + + +# ╔═╡ 05b63595-61dd-416f-8610-0f6ee1919a6f +@mytime res = Optim.optimize(f, g!, rec0, LBFGS(), + Optim.Options(iterations = 80, + store_trace=true)) + +# ╔═╡ 9ed4c43f-1ad5-43b1-bef6-40faf027207d +plot([t.iteration for t in res.trace], [t.value for t in res.trace], label="LBFGS", yaxis=:log) + +# ╔═╡ 5372b296-75d6-497d-b6c9-e2f104e76522 +[simshow(Array(res.minimizer) .* cispi(1.5), γ=1) simshow(Array(Array(beam)), γ=1)] + +# ╔═╡ 82273df3-554e-43f9-a881-63609c239648 +cat(simshow(Array(res.minimizer) .* cispi(1.5), γ=1), simshow(Array(Array(beam)), γ=1), dims=3)[:, :, 1] + +# ╔═╡ 02feb5a2-d43d-4191-b9cd-8d3946c05332 +simshow(Array(beam)) + +# ╔═╡ 980e4c15-7a04-4aad-89f1-119c33693288 + + +# ╔═╡ Cell order: +# ╠═d551f366-0213-11ef-176b-7fcf5b47315e +# ╠═9e6d909f-06de-46e7-929a-61bfd67bea26 +# ╠═20c5faa9-5694-42fd-b1bb-0c5093f6c93c +# ╠═b8cf117e-d084-4d42-8022-03fd6ac4c417 +# ╠═ef79854f-0ae9-41aa-8f92-aae8d7023d71 +# ╠═f9a13d72-f904-467c-8346-a9ef272693ff +# ╠═94ea2891-db37-477e-b535-1a6d402d4af2 +# ╠═a971769f-1215-4422-a209-fad710be16dd +# ╠═bdf3b801-28e0-41cf-9dc4-92fe8fdaca0e +# ╠═36b798c9-feee-43b7-b9ec-4f0558a44d2c +# ╠═ad4ef631-841d-412e-a09a-4db60fa7e4a9 +# ╠═15011619-f26a-4858-b0a0-d5dbcad49824 +# ╠═48048a82-5bbf-4f16-a7d9-15993c83f905 +# ╠═d058c9ba-82ec-4f1e-8938-eb9e1e95b713 +# ╠═5d21e0bc-ea6e-4b2c-8977-d1761d1b3c65 +# ╟─3f791571-2a36-4479-a853-c6ed4e8a270d +# ╠═e52cfcd8-eec8-4a8b-b953-f1c6ef2be959 +# ╠═df72f725-8464-46b0-a291-8503ea4e563f +# ╠═42c00ce2-d1e4-4f84-84be-a7e4cd6551ce +# ╠═90d78713-0d4d-48b2-b590-a79d4e7c12d5 +# ╠═f9c4d5ca-a35a-462a-b5f4-4120012f5a1b +# ╠═d84457d3-30ae-4711-b47f-7d2a7376a237 +# ╠═d219be8b-2d22-4fbe-b022-217edfeec5dc +# ╠═692e0263-8273-45a6-b3e6-c30bbeaeb649 +# ╠═5c550775-a614-44dc-a81d-638448d23067 +# ╠═56a03baf-3873-43d5-b8c1-8569de6026b0 +# ╠═fc9bb298-5b46-4585-9ea5-4df59405d18b +# ╠═de16936c-7774-4aab-9702-f0cf4fa660fd +# ╠═ee5260f3-d3a4-4bf2-8b54-31723c3c516f +# ╠═d90bff7e-150d-407d-a510-385fd5a32df3 +# ╠═eadee3cb-1302-4f1a-a21c-17265d9a2e6b +# ╠═05b63595-61dd-416f-8610-0f6ee1919a6f +# ╠═9ed4c43f-1ad5-43b1-bef6-40faf027207d +# ╠═5372b296-75d6-497d-b6c9-e2f104e76522 +# ╠═82273df3-554e-43f9-a881-63609c239648 +# ╠═02feb5a2-d43d-4191-b9cd-8d3946c05332 +# ╠═980e4c15-7a04-4aad-89f1-119c33693288 diff --git a/examples/Manifest.toml b/examples/Manifest.toml index 667affc..4993bc5 100644 --- a/examples/Manifest.toml +++ b/examples/Manifest.toml @@ -54,14 +54,15 @@ version = "0.4.0" [[deps.ArrayInterface]] deps = ["Adapt", "LinearAlgebra", "SparseArrays", "SuiteSparse"] -git-tree-sha1 = "44691067188f6bd1b2289552a23e4b7572f4528d" +git-tree-sha1 = "133a240faec6e074e07c31ee75619c90544179cf" uuid = "4fba245c-0d91-5ea0-9b3e-6abc04ee57a9" -version = "7.9.0" +version = "7.10.0" [deps.ArrayInterface.extensions] ArrayInterfaceBandedMatricesExt = "BandedMatrices" ArrayInterfaceBlockBandedMatricesExt = "BlockBandedMatrices" ArrayInterfaceCUDAExt = "CUDA" + ArrayInterfaceCUDSSExt = "CUDSS" ArrayInterfaceChainRulesExt = "ChainRules" ArrayInterfaceGPUArraysCoreExt = "GPUArraysCore" ArrayInterfaceReverseDiffExt = "ReverseDiff" @@ -72,6 +73,7 @@ version = "7.9.0" BandedMatrices = "aae01518-5342-5314-be14-df237901396f" BlockBandedMatrices = "ffab5731-97b5-5995-9138-79e8c1846df0" CUDA = "052768ef-5323-5732-b1bb-66c8b64840ba" + CUDSS = "45b445bb-4962-46a0-9369-b4df9d0f772e" ChainRules = "082447d4-558c-5d27-93f4-14fc19e9eca2" GPUArraysCore = "46192b85-c4d5-4398-a991-12ede77f4527" ReverseDiff = "37e2e3b7-166d-5795-8a7a-e32c996b4267" @@ -157,9 +159,9 @@ version = "0.5.0" [[deps.CUDA]] deps = ["AbstractFFTs", "Adapt", "BFloat16s", "CEnum", "CUDA_Driver_jll", "CUDA_Runtime_Discovery", "CUDA_Runtime_jll", "Crayons", "DataFrames", "ExprTools", "GPUArrays", "GPUCompiler", "KernelAbstractions", "LLVM", "LLVMLoopInfo", "LazyArtifacts", "Libdl", "LinearAlgebra", "Logging", "NVTX", "Preferences", "PrettyTables", "Printf", "Random", "Random123", "RandomNumbers", "Reexport", "Requires", "SparseArrays", "StaticArrays", "Statistics"] -git-tree-sha1 = "3dcab8a2c18ca319ea15a41d90e9528b8e93894a" +git-tree-sha1 = "dd1c682b372b6791b69f6823afe364fc92a0146c" uuid = "052768ef-5323-5732-b1bb-66c8b64840ba" -version = "5.3.0" +version = "5.3.1" weakdeps = ["ChainRulesCore", "SpecialFunctions"] [deps.CUDA.extensions] @@ -550,9 +552,9 @@ uuid = "7b1f6079-737a-58dc-b8bc-7a2ca5c1b5ee" [[deps.FillArrays]] deps = ["LinearAlgebra"] -git-tree-sha1 = "bfe82a708416cf00b73a3198db0859c82f741558" +git-tree-sha1 = "881275fc6b8c6f0dfb9cfa4a878979a33cb26be3" uuid = "1a297f60-69ca-5386-bcde-b61e274b549b" -version = "1.10.0" +version = "1.10.1" [deps.FillArrays.extensions] FillArraysPDMatsExt = "PDMats" @@ -726,9 +728,9 @@ version = "1.0.2" [[deps.HTTP]] deps = ["Base64", "CodecZlib", "ConcurrentUtilities", "Dates", "ExceptionUnwrapping", "Logging", "LoggingExtras", "MbedTLS", "NetworkOptions", "OpenSSL", "Random", "SimpleBufferStream", "Sockets", "URIs", "UUIDs"] -git-tree-sha1 = "8e59b47b9dc525b70550ca082ce85bcd7f5477cd" +git-tree-sha1 = "2c3ec1f90bb4a8f7beafb0cffea8a4c3f4e636ab" uuid = "cd3eb016-35fb-5094-929b-558a96fad6f3" -version = "1.10.5" +version = "1.10.6" [[deps.HarfBuzz_jll]] deps = ["Artifacts", "Cairo_jll", "Fontconfig_jll", "FreeType2_jll", "Glib_jll", "Graphite2_jll", "JLLWrappers", "Libdl", "Libffi_jll", "Pkg"] @@ -756,9 +758,9 @@ version = "0.2.4" [[deps.IRTools]] deps = ["InteractiveUtils", "MacroTools", "Test"] -git-tree-sha1 = "5d8c5713f38f7bc029e26627b687710ba406d0dd" +git-tree-sha1 = "d05027a62b4c9a2223820a9fdeae1110ad3946a5" uuid = "7869d1d1-7146-5819-86e3-90919afe41df" -version = "0.4.12" +version = "0.4.13" [[deps.IfElse]] git-tree-sha1 = "debdd00ffef04665ccbb3e150747a77560e8fad1" @@ -902,9 +904,9 @@ version = "0.1.5" [[deps.IntelOpenMP_jll]] deps = ["Artifacts", "JLLWrappers", "Libdl"] -git-tree-sha1 = "5fdf2fe6724d8caabf43b557b84ce53f3b7e2f6b" +git-tree-sha1 = "be50fe8df3acbffa0274a744f1a99d29c45a57f4" uuid = "1d5cc7b8-4909-519e-a0f8-d0f5ad9712d0" -version = "2024.0.2+0" +version = "2024.1.0+0" [[deps.InteractiveUtils]] deps = ["Markdown"] @@ -1219,10 +1221,10 @@ uuid = "6c6e2e6c-3030-632d-7369-2d6c69616d65" version = "0.1.4" [[deps.MKL_jll]] -deps = ["Artifacts", "IntelOpenMP_jll", "JLLWrappers", "LazyArtifacts", "Libdl"] -git-tree-sha1 = "72dc3cf284559eb8f53aa593fe62cb33f83ed0c0" +deps = ["Artifacts", "IntelOpenMP_jll", "JLLWrappers", "LazyArtifacts", "Libdl", "oneTBB_jll"] +git-tree-sha1 = "80b2833b56d466b3858d565adcd16a4a05f2089b" uuid = "856f044c-d86e-5d09-b602-aeab76dc8ba7" -version = "2024.0.0+0" +version = "2024.1.0+0" [[deps.MLStyle]] git-tree-sha1 = "bc38dff0548128765760c79eb7388a4b37fae2c8" @@ -1329,9 +1331,9 @@ version = "7.8.3" [[deps.NNlib]] deps = ["Adapt", "Atomix", "ChainRulesCore", "GPUArraysCore", "KernelAbstractions", "LinearAlgebra", "Pkg", "Random", "Requires", "Statistics"] -git-tree-sha1 = "1fa1a14766c60e66ab22e242d45c1857c83a3805" +git-tree-sha1 = "5055845dd316575ae2fc1f6dcb3545ff15fe547a" uuid = "872c559c-99b0-510c-b3b7-b6c96a88d5cd" -version = "0.9.13" +version = "0.9.14" [deps.NNlib.extensions] NNlibAMDGPUExt = "AMDGPU" @@ -1448,9 +1450,9 @@ version = "0.8.1+2" [[deps.OpenSSL]] deps = ["BitFlags", "Dates", "MozillaCACerts_jll", "OpenSSL_jll", "Sockets"] -git-tree-sha1 = "af81a32750ebc831ee28bdaaba6e1067decef51e" +git-tree-sha1 = "38cb508d080d21dc1128f7fb04f20387ed4c0af4" uuid = "4d8831e6-92b7-49fb-bdf8-b643e874388c" -version = "1.4.2" +version = "1.4.3" [[deps.OpenSSL_jll]] deps = ["Artifacts", "JLLWrappers", "Libdl"] @@ -1590,9 +1592,9 @@ version = "1.0.4" [[deps.PlutoUI]] deps = ["AbstractPlutoDingetjes", "Base64", "ColorTypes", "Dates", "FixedPointNumbers", "Hyperscript", "HypertextLiteral", "IOCapture", "InteractiveUtils", "JSON", "Logging", "MIMEs", "Markdown", "Random", "Reexport", "URIs", "UUIDs"] -git-tree-sha1 = "71a22244e352aa8c5f0f2adde4150f62368a3f2e" +git-tree-sha1 = "ab55ee1510ad2af0ff674dbcced5e94921f867a9" uuid = "7f904dfe-b85e-4ff6-b463-dae2292396a8" -version = "0.7.58" +version = "0.7.59" [[deps.PooledArrays]] deps = ["DataAPI", "Future"] @@ -2434,6 +2436,12 @@ deps = ["Artifacts", "Libdl"] uuid = "8e850ede-7688-5339-a07c-302acd2aaf8d" version = "1.52.0+1" +[[deps.oneTBB_jll]] +deps = ["Artifacts", "JLLWrappers", "Libdl"] +git-tree-sha1 = "7d0ea0f4895ef2f5cb83645fa689e52cb55cf493" +uuid = "1317d2d5-d96f-522e-a858-c73665f53c3e" +version = "2021.12.0+0" + [[deps.p7zip_jll]] deps = ["Artifacts", "Libdl"] uuid = "3f19e933-33d8-53b3-aaab-bd5110c3b7a0"