Skip to content

Commit 4773825

Browse files
a2raya2ray
authored
Merge pull request #74 from GeoscienceAustralia/silentlog
Silentlog
2 parents f3063ff + 9954a09 commit 4773825

17 files changed

+254
-150
lines changed

‎Project.toml

Lines changed: 3 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -1,7 +1,7 @@
11
name = "HiQGA"
22
uuid = "01574e87-63b6-4466-925a-334cf7e21b6b"
33
authors = ["richardt94 <29562583+richardt94@users.noreply.github.com>"]
4-
version = "0.4.7"
4+
version = "0.4.8"
55

66
[deps]
77
CSV = "336ed68f-0bac-5ca0-87d4-7b16caf5d00b"
@@ -25,6 +25,7 @@ KernelDensitySJ = "49dc5b4e-6806-4504-b2cc-a81764e7a38b"
2525
LazyGrids = "7031d0ef-c40d-4431-b2f8-61a8d2f650db"
2626
LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
2727
LinearMaps = "7a12625a-238d-50fd-b39a-03d52299707e"
28+
Logging = "56ddb016-857b-54e1-b83d-db4d58db5568"
2829
NearestNeighbors = "b8a86587-4115-5ab1-83bc-aa920d37bbce"
2930
Optim = "429524aa-4258-5aef-a3af-852621145aeb"
3031
Pkg = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f"
@@ -50,7 +51,7 @@ Dates = "1.7"
5051
DelimitedFiles = "1.7"
5152
Distances = "0.10.7"
5253
Distributed = "1.7"
53-
DistributedArrays = "0.6.6"
54+
DistributedArrays = "0.6.7"
5455
Distributions = "0.25.58"
5556
FastGaussQuadrature = "0.4.9"
5657
FileIO = "1.10"

‎examples/VTEM/DarlingParoo/McMC/04_plot.jl

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,4 @@
1-
dr = 20
1+
dr = 10
22
idx = [1,2]
33
burninfrac = 0.5
44
vmin, vmax = -2.5, 0.5

‎examples/tempest/Menindee/McMC/01_read_data.jl

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -25,5 +25,5 @@ soundings = transD_GP.TEMPEST1DInversion.read_survey_files(
2525
startfrom = 1,
2626
skipevery = 1,
2727
linenum = 1,
28-
dotillsounding = 2,
28+
dotillsounding = 4,
2929
makeqcplots = true)

‎examples/tempest/Menindee/McMC/02_set_options.jl

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,5 @@
11
## same for all soundingss
2-
nsamples = 201
2+
nsamples = 3001
33
nchainspersounding = 3
44
ppn = 4
55
@assert mod(ppn, nchainspersounding+1) == 0

‎examples/tempest/Menindee/McMC/04_plot.jl

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -3,7 +3,7 @@ idxplot = [1,2]
33
burninfrac = 0.5
44
nbins = 50
55
computeforwards, nforwards = true, 2
6-
interpolatedr = 10
6+
interpolatedr = 1
77
vmax, vmin = .- extrema(fbounds)
88
## get the summary images of resistivity and nuisances over the profile
99
transD_GP.AEMwithNuisanceMcMCInversionTools.summaryAEMwithnuisanceimages(soundings, opt, optn; zall,

‎src/AEM_VMD_HMD.jl

Lines changed: 4 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -233,8 +233,10 @@ function checkrampformintime(times, ramp, minresptime, maxtime, rampchoice, isdI
233233
end
234234
end
235235
end
236-
mintime = max(0.5minta, minresptime)
237-
# though I believe mintime = minta always is safe.
236+
mintime = max(0.25minta, minresptime)
237+
# I believe mintime = 0.25minta is safe.
238+
# later times can cause errors in tb, unless extrapolating earlier
239+
# than F.interptimes[1]. Extrapolation I'm hesitant to do.
238240
# if maxta > maxtime
239241
maxtime = min(1.5maxta, maxtime)
240242
# 1.5 to make sure we clear maxta in the interptimes

‎src/AEMnoNuisanceGradientInversionTools.jl

Lines changed: 24 additions & 22 deletions
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,6 @@
11
module AEMnoNuisanceGradientInversionTools
22
using Distributed, Dates, Printf, PyPlot, DelimitedFiles, StatsBase
3-
using ..AbstractOperator, ..CommonToAll, ..GP
3+
using ..AbstractOperator, ..CommonToAll, ..GP, ..SoundingDistributor
44
import ..AbstractOperator.makeoperator
55
import ..AbstractOperator.Sounding
66
import ..AbstractOperator.returnforwrite
@@ -15,26 +15,26 @@ function compress(soundings, zall; prefix="", rmfile=true, isfirstparalleliterat
1515
fname = soundings[1].sounding_string*"_gradientinv.dat"
1616
!isfile(fname) && throw(AssertionError("file does not exist perhaps soundings already zipped?"))
1717
fout = prefix == "" ? "zipped.dat" : prefix*"_zipped.dat"
18-
iomode = "w"
19-
if isfile(fout)
20-
isfirstparalleliteration && throw(AssertionError("Zipped file "*fout*" exists, will not overwrite!"))
21-
iomode = "a"
22-
end
23-
io = open(fout, iomode)
18+
isfirstparalleliteration && isfile(fout) && throw(AssertionError("Zipped file "*fout*" exists, will not overwrite!"))
2419
for (i, s) in enumerate(soundings)
2520
fname = s.sounding_string*"_gradientinv.dat"
2621
A = readdlm(fname)
2722
ϕd = A[end,2]
2823
σgrid = vec(A[end,3:end])
29-
for el in [returnforwrite(s)...; vec(zall); σgrid; ϕd]
30-
msg = @sprintf("%.4f\t", el)
31-
write(io, msg)
32-
end
33-
write(io, "\n")
34-
flush(io) # slower but ensures write is complete
24+
elinonerow = [returnforwrite(s)..., vec(zall), σgrid, ϕd]
25+
nelinonerow = length(elinonerow)
26+
writenames = [string.(s.writefields)..., "zcenter", "log10_cond", "phid_err"]
27+
sfmt = fill("%15.3f", nelinonerow)
28+
ϕd > 1e4 && (ϕd = 1e4 )
29+
iomode = (isfirstparalleliteration && i==1) ? "w" : "a"
30+
writeasegdat(elinonerow, sfmt, fout[1:end-4], iomode)
3531
rmfile && rm(fname)
32+
if isfirstparalleliteration && i == 1
33+
channel_names = [writenames, fill("", nelinonerow), writenames]
34+
writeasegdfnfromonerow(elinonerow, channel_names, sfmt, fout[1:end-4])
35+
dfn2hdr(fout[1:end-4]*".dfn")
36+
end
3637
end
37-
close(io)
3838
end
3939

4040
# plot the convergence and the result
@@ -237,21 +237,22 @@ function loopacrossAEMsoundings(soundings::Array{S, 1}, aem_in, σstart, σ0;
237237
compresssoundings = true,
238238
zipsaveprefix = "",
239239
minimprovfrac = nothing,
240+
verbose = false,
240241
minimprovkickinstep = round(Int, nstepsmax/2),
241242
) where S<:Sounding
242243

243244
@assert nsequentialiters != -1
244245
nparallelsoundings = nworkers()
245246
nsoundings = length(soundings)
246247
zall = zboundarytocenter(aem_in.z[aem_in.nfixed+1:end]) # needed for sounding compression in write
248+
249+
writetogloballog("will require $nsequentialiters iterations of $nparallelsoundings soundings", iomode="w")
250+
writetogloballog("starting sequential parallel iterations at $(Dates.now())")
247251
for iter = 1:nsequentialiters
248-
if iter<nsequentialiters
249-
ss = (iter-1)*nparallelsoundings+1:iter*nparallelsoundings
250-
else
251-
ss = (iter-1)*nparallelsoundings+1:nsoundings
252-
end
253-
@info "soundings in loop $iter of $nsequentialiters", ss
252+
ss = getss_deterministic(iter, nsequentialiters, nparallelsoundings, nsoundings)
253+
writetogloballog("soundings in loop $iter of $nsequentialiters $ss")
254254
pids = workers()
255+
t2 = time()
255256
@sync for (i, s) in enumerate(ss)
256257
aem = makeoperator(aem_in, soundings[s])
257258
fname = soundings[s].sounding_string*"_gradientinv.dat"
@@ -279,15 +280,16 @@ function loopacrossAEMsoundings(soundings::Array{S, 1}, aem_in, σstart, σ0;
279280
breakonknown = breakonknown ,
280281
dobo = dobo ,
281282
fname = fname ,
282-
minimprovfrac,
283+
minimprovfrac, verbose,
283284
minimprovkickinstep)
284285

285286

286287
end # @sync
287288
isfirstparalleliteration = iter == 1 ? true : false
288289
compresssoundings && compress(soundings[ss[1]:ss[end]], zall,
289290
isfirstparalleliteration = isfirstparalleliteration, prefix=zipsaveprefix)
290-
@info "done $iter out of $nsequentialiters at $(Dates.now())"
291+
dt = time() - t2 # seconds
292+
writetogloballog("done $iter out of $nsequentialiters at $(Dates.now()) in $dt sec")
291293
end
292294
end
293295

‎src/AEMnoNuisanceMcMCInversionTools.jl

Lines changed: 11 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -37,6 +37,7 @@ function summaryAEMimages(soundings::Array{S, 1}, opt::Options;
3737
yl = nothing,
3838
showplot = true,
3939
showmean = false,
40+
Rmax = nothing,
4041
logscale = true,
4142
dpi = 300) where S<:Sounding
4243
compatidxwarn(idx, lnames)
@@ -48,7 +49,7 @@ function summaryAEMimages(soundings::Array{S, 1}, opt::Options;
4849
continueflag && continue
4950
summaryimages(soundings[a:b], opt; qp1, qp2, burninfrac, zall,dz, dr,
5051
fontsize, vmin, vmax, cmap, figsize, topowidth, idx=idspec, omitconvergence, useML,
51-
preferEright, showplot, preferNright, saveplot, yl, dpi, showmean, logscale)
52+
preferEright, showplot, preferNright, saveplot, yl, dpi, showmean, logscale, Rmax)
5253
end
5354
nothing
5455
end
@@ -76,7 +77,8 @@ function summaryimages(soundings::Array{S, 1}, opt::Options;
7677
logscale = true,
7778
showplot = true,
7879
showmean = false,
79-
dpi = 300) where S<:Sounding
80+
dpi = 300,
81+
Rmax=nothing) where S<:Sounding
8082
@assert !(preferNright && preferEright) # can't prefer both labels to the right
8183
pl, pm, ph, ρmean, χ²mean, χ²sd = summarypost(soundings, opt; zall, qp1, qp2, burninfrac, useML)
8284

@@ -85,9 +87,10 @@ function summaryimages(soundings::Array{S, 1}, opt::Options;
8587
zall, dz; dr)
8688

8789
lname = "Line_$(soundings[1].linenum)"
90+
8891
plotsummarygrids1(soundings, σmeangrid, phgrid, plgrid, pmgrid, gridx, gridz, topofine, R, Z, χ²mean, χ²sd, lname; qp1, qp2,
8992
figsize, fontsize, cmap, vmin, vmax,
90-
topowidth, idx, omitconvergence, useML,
93+
topowidth, idx, omitconvergence, useML, Rmax,
9194
preferEright, preferNright, saveplot, showplot, dpi,
9295
yl, showmean, logscale)
9396
end
@@ -145,7 +148,7 @@ function processonesounding(opt_in::Options, sounding::Sounding, zall, burninfra
145148
ndata = getndata(sounding)
146149
χ²mean = mean(χ²)/ndata
147150
χ²sd = std(χ²)/ndata
148-
if hasproperty(sounding, :force_ML)
151+
if hasproperty(sounding, :forceML)
149152
if sounding.forceML
150153
χ²mean = 1.
151154
χ²sd = 0.
@@ -334,9 +337,10 @@ function loopacrossAEMsoundings(soundings::Array{S, 1}, aem_in::Operator1D, opt_
334337
nsoundings = length(soundings)
335338
nsequentialiters, nparallelsoundings = splittasks(soundings; nchainspersounding, ppn)
336339

340+
writetogloballog("starting sequential parallel iterations at $(Dates.now())")
337341
for iter = 1:nsequentialiters
338342
ss = getss(iter, nsequentialiters, nparallelsoundings, nsoundings)
339-
@info "soundings in loop $iter of $nsequentialiters", ss
343+
writetogloballog("soundings in loop $iter of $nsequentialiters $ss")
340344
t2 = time()
341345
@sync for (i, s) in Iterators.reverse(enumerate(ss))
342346
pids = getpids(i, nchainspersounding)
@@ -350,8 +354,8 @@ function loopacrossAEMsoundings(soundings::Array{S, 1}, aem_in::Operator1D, opt_
350354

351355
end # @sync
352356
dt = time() - t2 #seconds
353-
t2 = time()
354-
@info "done $iter out of $nsequentialiters at $(Dates.now()) in $dt sec"
357+
catlocallogs(nparallelsoundings, nchainspersounding)
358+
writetogloballog("done $iter out of $nsequentialiters at $(Dates.now()) in $dt sec")
355359
end
356360
end
357361

‎src/AEMwithNuisanceGradientInversionTools.jl

Lines changed: 27 additions & 22 deletions
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,6 @@
11
module AEMwithNuisanceGradientInversionTools
22
using Distributed, Dates, Printf, PyPlot, DelimitedFiles, StatsBase
3-
using ..AbstractOperator, ..CommonToAll, ..GP
3+
using ..AbstractOperator, ..CommonToAll, ..GP, ..SoundingDistributor
44
import ..AbstractOperator.makeoperator
55
import ..AbstractOperator.setnuboundsandstartforgradinv
66
import ..AbstractOperator.Sounding
@@ -17,27 +17,27 @@ function compress(soundings, zall, nnu; prefix="", rmfile=true, isfirstparalleli
1717
fname = soundings[1].sounding_string*"_gradientinv.dat"
1818
!isfile(fname) && throw(AssertionError("file does not exist perhaps soundings already zipped?"))
1919
fout = prefix == "" ? "zipped.dat" : prefix*"_zipped.dat"
20-
iomode = "w"
21-
if isfile(fout)
22-
isfirstparalleliteration && throw(AssertionError("Zipped file "*fout*" exists, will not overwrite!"))
23-
iomode = "a"
24-
end
25-
io = open(fout, iomode)
20+
isfirstparalleliteration && isfile(fout) && throw(AssertionError("Zipped file "*fout*" exists, will not overwrite!"))
2621
for (i, s) in enumerate(soundings)
2722
fname = s.sounding_string*"_gradientinv.dat"
2823
A = readdlm(fname)
2924
ϕd = A[end,2]
3025
σgrid = vec(A[end,3:end-nnu])
3126
nu = vec(A[end,end-nnu+1:end])
32-
for el in [returnforwrite(s)...; vec(zall); σgrid; nu; ϕd]
33-
msg = @sprintf("%.4f\t", el)
34-
write(io, msg)
35-
end
36-
write(io, "\n")
37-
flush(io) # slower but ensures write is complete
27+
elinonerow = [returnforwrite(s)..., vec(zall), σgrid, nu, ϕd]
28+
nelinonerow = length(elinonerow)
29+
writenames = [string.(s.writefields)..., "zcenter", "log10_cond", "nuisance_param_in_order", "phid_err"]
30+
sfmt = fill("%15.3f", nelinonerow)
31+
ϕd > 1e4 && (ϕd = 1e4 )
32+
iomode = (isfirstparalleliteration && i==1) ? "w" : "a"
33+
writeasegdat(elinonerow, sfmt, fout[1:end-4], iomode)
3834
rmfile && rm(fname)
35+
if isfirstparalleliteration && i == 1
36+
channel_names = [writenames, fill("", nelinonerow), writenames]
37+
writeasegdfnfromonerow(elinonerow, channel_names, sfmt, fout[1:end-4])
38+
dfn2hdr(fout[1:end-4]*".dfn")
39+
end
3940
end
40-
close(io)
4141
end
4242

4343
# plot the convergence and the result
@@ -234,6 +234,9 @@ function loopacrossAEMsoundings(soundings::Array{S, 1}, aem_in, σstart, σ0, nu
234234
knownvalue = 0.7,
235235
compresssoundings = true,
236236
zipsaveprefix = "",
237+
minimprovfrac = nothing,
238+
verbose = false,
239+
minimprovkickinstep = round(Int, nstepsmax/2),
237240
# optim stuff
238241
ntriesnu = 5,
239242
boxiters = 3,
@@ -247,14 +250,14 @@ function loopacrossAEMsoundings(soundings::Array{S, 1}, aem_in, σstart, σ0, nu
247250
nparallelsoundings = nworkers()
248251
nsoundings = length(soundings)
249252
zall = zboundarytocenter(aem_in.z[aem_in.nfixed+1:end]) # needed for sounding compression in write
253+
254+
writetogloballog("will require $nsequentialiters iterations of $nparallelsoundings soundings", iomode="w")
255+
writetogloballog("starting sequential parallel iterations at $(Dates.now())")
250256
for iter = 1:nsequentialiters
251-
if iter<nsequentialiters
252-
ss = (iter-1)*nparallelsoundings+1:iter*nparallelsoundings
253-
else
254-
ss = (iter-1)*nparallelsoundings+1:nsoundings
255-
end
256-
@info "soundings in loop $iter of $nsequentialiters", ss
257+
ss = getss_deterministic(iter, nsequentialiters, nparallelsoundings, nsoundings)
258+
writetogloballog("soundings in loop $iter of $nsequentialiters $ss")
257259
pids = workers()
260+
t2 = time()
258261
@sync for (i, s) in enumerate(ss)
259262
aem = makeoperator(aem_in, soundings[s])
260263
nubounds, nustart = setnuboundsandstartforgradinv(aem, nuboundsΔ)
@@ -273,14 +276,16 @@ function loopacrossAEMsoundings(soundings::Array{S, 1}, aem_in, σstart, σ0, nu
273276
regularizeupdate = regularizeupdate,
274277
knownvalue = knownvalue ,
275278
fname = fname ,
276-
ntriesnu, nubounds, boxiters, usebox, reducenuto, debuglevel, breaknuonknown)
279+
ntriesnu, nubounds, boxiters, usebox, reducenuto, debuglevel, breaknuonknown,
280+
minimprovfrac, verbose, minimprovkickinstep)
277281

278282

279283
end # @sync
280284
isfirstparalleliteration = iter == 1 ? true : false
281285
compresssoundings && compress(soundings[ss[1]:ss[end]], zall, size(nuboundsΔ, 1),
282286
isfirstparalleliteration = isfirstparalleliteration, prefix=zipsaveprefix)
283-
@info "done $iter out of $nsequentialiters at $(Dates.now())"
287+
dt = time() - t2 # seconds
288+
writetogloballog("done $iter out of $nsequentialiters at $(Dates.now()) in $dt sec")
284289
end
285290
end
286291

‎src/AEMwithNuisanceMcMCInversionTools.jl

Lines changed: 4 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -173,9 +173,10 @@ function loopacrossAEMsoundings(soundings::Array{S, 1}, aem_in::Operator1D, opt_
173173
nsoundings = length(soundings)
174174
nsequentialiters, nparallelsoundings = splittasks(soundings; nchainspersounding, ppn)
175175

176+
writetogloballog("starting sequential parallel iterations at $(Dates.now())")
176177
for iter = 1:nsequentialiters
177178
ss = getss(iter, nsequentialiters, nparallelsoundings, nsoundings)
178-
@info "soundings in loop $iter of $nsequentialiters", ss
179+
writetogloballog("soundings in loop $iter of $nsequentialiters $ss")
179180
t2 = time()
180181
@sync for (i, s) in Iterators.reverse(enumerate(ss))
181182
pids = getpids(i, nchainspersounding)
@@ -190,8 +191,8 @@ function loopacrossAEMsoundings(soundings::Array{S, 1}, aem_in::Operator1D, opt_
190191

191192
end # @sync
192193
dt = time() - t2 #seconds
193-
t2 = time()
194-
@info "done $iter out of $nsequentialiters at $(Dates.now()) in $dt sec"
194+
catlocallogs(nparallelsoundings, nchainspersounding)
195+
writetogloballog("done $iter out of $nsequentialiters at $(Dates.now()) in $dt sec")
195196
end
196197
end
197198

0 commit comments

Comments
 (0)