Skip to content

Commit 16bdb6e

Browse files
committed
Add tests and docstrings
1 parent d542c64 commit 16bdb6e

File tree

6 files changed

+191
-27
lines changed

6 files changed

+191
-27
lines changed

docs/src/api.md

Lines changed: 2 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -6,10 +6,6 @@
66
Pages = ["api.md"]
77
```
88

9-
```@docs
10-
Comrade.Comrade
11-
```
12-
139
## Model Definitions
1410

1511

@@ -133,13 +129,13 @@ Comrade.dataproducts
133129
Comrade.skymodel
134130
Comrade.instrumentmodel(::Comrade.AbstractVLBIPosterior)
135131
Comrade.instrumentmodel(::Comrade.AbstractVLBIPosterior, ::Any)
136-
Comrade.forward_model
137132
Comrade.prior_sample
138133
Comrade.likelihood
139134
Comrade.VLBIPosterior
140135
Comrade.simulate_observation
141136
Comrade.residuals
142137
Comrade.TransformedVLBIPosterior
138+
Comrade.ImgNormalData
143139
HypercubeTransform.transform(::Comrade.TransformedVLBIPosterior, ::Any)
144140
HypercubeTransform.inverse(::Comrade.TransformedVLBIPosterior, ::Any)
145141
HypercubeTransform.ascube(::Comrade.VLBIPosterior)
@@ -181,6 +177,7 @@ Comrade.rmap
181177
```@docs
182178
Comrade.build_datum
183179
Comrade.ObservedSkyModel
180+
Comrade.forward_model_map
184181
```
185182

186183
### eht-imaging interface (Internal)

docs/src/base_api.md

Lines changed: 14 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -62,9 +62,10 @@ ComradeBase.intensitymap_numeric!
6262

6363
### Image Domain
6464
```@docs
65-
ComradeBase.imagepixels
65+
ComradeBase.AbstractDualDomain
6666
ComradeBase.RectiGrid
6767
ComradeBase.UnstructuredDomain
68+
ComradeBase.imagepixels
6869
ComradeBase.dims
6970
ComradeBase.named_dims
7071
ComradeBase.axisdims
@@ -89,8 +90,20 @@ ComradeBase.baseimage
8990
ComradeBase.centroid
9091
ComradeBase.second_moment
9192
ComradeBase.stokes
93+
ComradeBase.dualmap
94+
ComradeBase.DualMap
95+
```
96+
97+
### Time and Frequency Domain
98+
```@docs
99+
ComradeBase.DomainParams
100+
ComradeBase.paramtype
101+
ComradeBase.getparam
102+
ComradeBase.@unpack_params
103+
ComradeBase.build_param
92104
```
93105

106+
94107
## Internal Methods not part of public API
95108
```@docs
96109
ComradeBase._visibilitymap

src/posterior/likelihood.jl

Lines changed: 20 additions & 19 deletions
Original file line numberDiff line numberDiff line change
@@ -106,43 +106,44 @@ function makelikelihood(data::Comrade.EHTObservationTable{<:Comrade.EHTClosurePh
106106
return
107107
end
108108

109+
export ImgNormalData
109110

110-
struct CentroidData{M, N}
111+
"""
112+
ImgNormalData(reduction, measurement, noise)
113+
114+
Container for image-domain data assuming a normal likelihood with a `measurement` and `noise`.
115+
"""
116+
struct ImgNormalData{F, M, N}
117+
reduction::F
111118
measurement::M
112119
noise::N
113120
end
114121

115122

116-
struct _Centroid{F, T}
123+
struct _Reduced{F, T}
117124
f::F
118-
Σ::T
125+
σ::T
119126
end
120127

121128
struct NormalFast{T, S}
122129
μ::T
123-
Σ::S
130+
σ::S
124131
end
125132

126133
function DensityInterface.logdensityof(d::NormalFast, x)
127-
dev = x .- d.μ
128-
exponent = - sum(abs2, dev) * inv(2 * d.Σ)
134+
dev = (x .- d.μ).*inv.(d.σ)
135+
exponent = -sum(abs2, dev)/2
129136
return exponent
130137
end
131138

132-
function (c::_Centroid)(μ)
133-
return NormalFast(c.f(μ), c.Σ)
139+
function (c::_Reduced)(μ)
140+
return NormalFast(c.f(μ), c.σ)
134141
end
135142

136-
function makelikelihood(data::CentroidData)
137-
Σ = data.noise .^ 2
143+
function makelikelihood(data::ImgNormalData)
144+
σ = data.noise
138145
meas = data.measurement
139-
f = SVector centroid
140-
= ConditionedLikelihood(_Centroid(f, Σ), SVector(meas))
146+
f = data.reduction
147+
= ConditionedLikelihood(_Reduced(f, σ), meas)
141148
return
142-
end
143-
144-
# function Comrade.centroid(m::VLBISkyModels.ModifiedModel{<:ContinuousImage}, g)
145-
# c = centroid(VLBISkyModels.unmodified(m), g)
146-
# tr = m.transforms
147-
# return shift_centroid(tr, c)
148-
# end
149+
end

test/Core/bayes.jl

Lines changed: 0 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -208,8 +208,6 @@ end
208208

209209
c2 = chi2(postsim, x; reduce = true)
210210
c2nn = chi2(postsim_nn, x; reduce = true)
211-
@info c2
212-
@info logdensityof(postsim, x)
213211
@test all(x -> reduce(&, x .< 1.25), c2)
214212
@test all(x -> reduce(&, x .≈ 0), c2nn)
215213

test/Core/core.jl

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -10,4 +10,5 @@ using VLBIImagePriors
1010
include(joinpath(@__DIR__, "observation.jl"))
1111
include(joinpath(@__DIR__, "partially_fixed.jl"))
1212
include(joinpath(@__DIR__, "models.jl"))
13+
include(joinpath(@__DIR__, "imgnormal.jl"))
1314
include(joinpath(@__DIR__, "bayes.jl"))

test/Core/imgnormal.jl

Lines changed: 154 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,154 @@
1+
using Distributions
2+
3+
@testset "ImgNormalData" begin
4+
# Load test data to get a working model
5+
_, vis, amp, lcamp, cphase = load_data()
6+
7+
# Create a simple image grid
8+
g = imagepixels(μas2rad(150.0), μas2rad(150.0), 32, 32)
9+
skym = SkyModel(test_model, test_prior(), g)
10+
11+
# Create a base posterior to sample from
12+
base_post = VLBIPosterior(skym, vis)
13+
14+
# Test with a simple reduction function (total flux)
15+
@testset "Total Flux Measurement" begin
16+
# Create reduction function for total flux
17+
reduction = flux
18+
19+
# Get a sample image to create mock measurement
20+
x0 = prior_sample(base_post)
21+
m0 = skymodel(base_post, x0)
22+
img0 = intensitymap(m0, g)
23+
24+
# Create mock measurement with noise
25+
true_flux = flux(img0)
26+
measurement = [true_flux]
27+
noise = [0.1]
28+
29+
# Create ImgNormalData object
30+
imgdata = ImgNormalData(reduction, measurement, noise)
31+
32+
# Test that fields are correctly stored
33+
@test imgdata.reduction === reduction
34+
@test imgdata.measurement == measurement
35+
@test imgdata.noise == noise
36+
37+
# Create likelihood from the image data
38+
= Comrade.makelikelihood(imgdata)
39+
40+
# Test that we can evaluate the likelihood
41+
@test ℓ isa Comrade.ConditionedLikelihood
42+
43+
# Test likelihood evaluation with the image
44+
lp = logdensityof(ℓ, img0)
45+
@test lp isa Real
46+
@test isfinite(lp)
47+
48+
# Test that likelihood is higher when image flux matches measurement
49+
# Create image with exact flux
50+
img_exact = img0 .* (measurement[1] / flux(img0))
51+
lp_exact = logdensityof(ℓ, img_exact)
52+
53+
# The exact match should have higher likelihood
54+
@test lp_exact >= lp
55+
end
56+
57+
@testset "Centroid Measurement" begin
58+
# Create reduction function for centroid
59+
reduction = img -> begin
60+
cx, cy = centroid(img)
61+
return [cx, cy]
62+
end
63+
64+
# Get a sample image
65+
x0 = prior_sample(base_post)
66+
m0 = skymodel(base_post, x0)
67+
img0 = intensitymap(m0, g)
68+
69+
# Create mock measurement
70+
cx0, cy0 = centroid(img0)
71+
measurement = [cx0, cy0]
72+
noise = [μas2rad(5.0), μas2rad(5.0)]
73+
74+
# Create ImgNormalData object
75+
imgdata = ImgNormalData(reduction, measurement, noise)
76+
77+
# Create likelihood
78+
= Comrade.makelikelihood(imgdata)
79+
80+
# Test likelihood evaluation
81+
lp = logdensityof(ℓ, img0)
82+
@test lp isa Real
83+
@test isfinite(lp)
84+
85+
# Test with a different image
86+
x1 = prior_sample(base_post)
87+
m1 = skymodel(base_post, x1)
88+
img1 = intensitymap(m1, g)
89+
lp1 = logdensityof(ℓ, img1)
90+
@test lp1 isa Real
91+
@test isfinite(lp1)
92+
end
93+
94+
@testset "Multiple Measurements" begin
95+
# Create reduction function that returns multiple quantities
96+
reduction = img -> begin
97+
f = flux(img)
98+
cx, cy = centroid(img)
99+
return [f, cx, cy]
100+
end
101+
102+
# Get a sample image
103+
x0 = prior_sample(base_post)
104+
m0 = skymodel(base_post, x0)
105+
img0 = intensitymap(m0, g)
106+
107+
# Create mock measurements
108+
f0 = flux(img0)
109+
cx0, cy0 = centroid(img0)
110+
measurement = [f0, cx0, cy0]
111+
noise = [0.1, μas2rad(5.0), μas2rad(5.0)]
112+
113+
# Create ImgNormalData object
114+
imgdata = ImgNormalData(reduction, measurement, noise)
115+
116+
# Create likelihood
117+
= Comrade.makelikelihood(imgdata)
118+
119+
# Test likelihood evaluation
120+
lp = logdensityof(ℓ, img0)
121+
@test lp isa Real
122+
@test isfinite(lp)
123+
124+
# Test that the likelihood is composed correctly
125+
# The reduction should return a vector
126+
red_result = reduction(img0)
127+
@test length(red_result) == 3
128+
@test red_result[1] f0
129+
@test red_result[2] cx0
130+
@test red_result[3] cy0
131+
end
132+
133+
@testset "Combining ImgNormalData with VLBI data" begin
134+
# Create image measurement data
135+
reduction = flux
136+
measurement = [1.0]
137+
noise = [0.1]
138+
imgdata = ImgNormalData(reduction, measurement, noise)
139+
140+
# Create posterior with both VLBI and image data using keyword argument
141+
post = VLBIPosterior(skym, vis; imgdata = (imgdata,))
142+
143+
# Test that we can evaluate the posterior
144+
x_prior = prior_sample(post)
145+
lp_post = logdensityof(post, x_prior)
146+
@test lp_post isa Real
147+
@test isfinite(lp_post)
148+
149+
# Test dataproducts
150+
dp = dataproducts(post)
151+
@test dp isa Tuple
152+
@test vis in dp
153+
end
154+
end

0 commit comments

Comments
 (0)