From fb7f811fa6db51e37899926e7ebd4d2bdfd5e81c Mon Sep 17 00:00:00 2001 From: Tom McGrath <72618709+tmcgrath325@users.noreply.github.com> Date: Tue, 4 Jun 2024 11:49:49 -0500 Subject: [PATCH] Add prepend! for ImagineSignal traces, stimulus control example (#91) * add prepend! functionality for editing traces * add workflow script for stimulus control * fix cases with negative stimulus delay * fix incomplete comment about injection delay * switch to pushfirst! * Fixed improper argument format for pushfirst! * Added tests for prepend! The tests are similar to the planned usage for prepend!, adding a constant value before an experiment * Changed example from bi- to unidirection --- examples/workflow_stim_control.jl | 132 ++++++++++++++++++++++++++++++ src/ImagineInterface.jl | 2 +- src/imaginesignal.jl | 43 ++++++++++ test/low_level.jl | 14 ++++ 4 files changed, 190 insertions(+), 1 deletion(-) create mode 100644 examples/workflow_stim_control.jl diff --git a/examples/workflow_stim_control.jl b/examples/workflow_stim_control.jl new file mode 100644 index 0000000..f1503b7 --- /dev/null +++ b/examples/workflow_stim_control.jl @@ -0,0 +1,132 @@ +using ImagineInterface + +import Unitful: μm, s + +# See "workflow.jl" for a more detailed overview of package functionality + +############PRIOR TO USE################# +# Note that before stimulus control via Imagine can be used, an appropriate Chromeleon program (.pgm file) must be used. +# See documentation elsewhere (Autosampler.jl) for instructions on creation of a .pgm file. + +############LOAD A RIG-SPECIFIC COMMAND TEMPLATE################# +sample_rate = 50000s^-1 +rig = "ocpi-2" +ocpi2 = rigtemplate(rig; sample_rate = sample_rate) + +pos = getpositioners(ocpi2)[1] #positioner trace +las1 = getlasers(ocpi2)[1] #laser trace +cam1 = getcameras(ocpi2)[1] #camera trace + +stim1 = getstimuli(ocpi2)[1] #stimulus timing trace + +############STACK PARAMETERS################# +pmin = 0.0*μm #Piezo start position in microns +pmax = 200.0*μm #Piezo stop position in microns +stack_img_time = 1.0s #Time to complete the imaging sweep with the piezo (remember we also need to reset it to its starting position) +reset_time = 0.5s #Only used for unidirectional sweeps. Time to reset piezo to starting position. This time plus "stack_img_time" determines how long it takes to complete an entire stack and be ready to start a new stack +z_spacing = 3.1μm #The space between slices in the z-stack. +z_pad = 5.0μm #Set this greater than 0 if you want to ignore the edges of the positioner sweep (only take slices in a central region) + +exp_time = 0.011s #Exposure time of the camera. Make sure this is greater than mn_exp and less than mx_exp above +flash_frac = 0.1 #fraction of time to keep laser on during exposure. If you set this greater than 1 then the laser will stay on constantly during the imaging sweep + +############STACK TIMING PARAMETERS################# +# These timing parameters depend on the specifics of your experiment. +# You have free control over the recording durations, but the stimulus lead time will depend on tubing volume between the autosampler and recording chamber + +stimulus_lead_time = 20s #Interval between injection time and stimulus arrival time (depends on flowrate and tubing volume, so you should measure this for your setup) +baseline_duration = 10s #Interval between recording start time and stimulus time (positive if recording starts before stimulus presentation) +stimulus_duration = 20s #Duration of recording for each trial (after stimulus presentation) + +# The total recording time for each trial is the sum of the baseline and stimulus durations +@show total_recording_duration = baseline_duration + stimulus_duration + +# We will also need to set the time between each trial, which should be enough to accommodate washing of the recording chamber and sample loop +inter_trial_duration = 30s #Duration of interval between trials + +# The total number of recordings to be obtained. This should be the product of the number of stimuli and number of replicates for your experiment +n_stimuli = 2 #Number of stimuli to present +n_replicates = 3 #Number of replicates for each stimulus +@show n_trials = n_stimuli * n_replicates #Number of trials to record + +############STACK GENERATION EXAMPLE################# +# unidirectional sweep waveforms +unidi_samps = gen_unidirectional_stack(pmin, pmax, z_spacing, stack_img_time, reset_time, exp_time, sample_rate, flash_frac; z_pad = z_pad) +append!(pos, "unidi_stack_pos", unidi_samps["positioner"]) +append!(las1, "unidi_stack_las1", unidi_samps["laser"]) +append!(cam1, "unidi_stack_cam1", unidi_samps["camera"]); +sweep_nframes = unidi_samps["nframes"]; #store this for later + +# Repeat the unidirectional waveform enough times to achieve the desired recording duration +@show sweeps_per_stim = Int(ceil(total_recording_duration / (stack_img_time + reset_time))) +replicate!(pos, (sweeps_per_stim-1)) +replicate!(las1, (sweeps_per_stim-1)) +replicate!(cam1, (sweeps_per_stim-1)) + +# inter-stimulus rest waveforms +wait_nsamps = inter_trial_duration * sample_rate +wait_samps = Dict( + "positioner" => fill(eltype(unidi_samps["positioner"])(pmin), wait_nsamps), + "laser" => fill(false, wait_nsamps), + "camera" => fill(false, wait_nsamps) +) +append!(pos, "wait_stack_pos", wait_samps["positioner"]) +append!(las1, "wait_stack_las1", wait_samps["laser"]) +append!(cam1, "wait_stack_cam1", wait_samps["camera"]); + +# We now have, for the positioner, camera, and laser traces, a waveform for recording a single trial + +# Injection trigger signal (true = on, false = off). The autosampler triggers an injection upon the switch to an "on" signal +stim_on_nsamps = Int(ceil(stimulus_duration * sample_rate)) #Number of samples for the "on" signal +stim_on_samps = fill(true, stim_on_nsamps) +stim_off_nsamps = Int(ceil((total_recording_duration + inter_trial_duration - stimulus_duration) * sample_rate)) +stim_off_samps = fill(false, stim_off_nsamps) + +stim1 = getstimuli(ocpi2)[1] +append!(stim1, "stim1", [stim_on_samps; stim_off_samps]) + +# We now have, for the stimulus trace, a waveform for triggering injection for a single trial + +############STIMULUS TIMING EXAMPLE################# +# Repeat each trace for the desired number of trials +replicate!(pos, n_trials-1) # Remember, we already have one trial's worth of data in the trace +replicate!(las1, n_trials-1) +replicate!(cam1, n_trials-1) +replicate!(stim1, n_trials-1) + +# Currently, the start of a recording and the injection trigger occur at the same time. +# We need to trigger the injection at the appropriate time to account for the desired baseline duration +# and the time it takes for the stimulus to arrive at the recording chamber for each trial +@show stim_delay = baseline_duration - stimulus_lead_time + +# Pad beginning and end of experiment to account for interval between recording start and injection trigger +pad_nsamps = Int(abs(stim_delay) * sample_rate) +pad_samps = Dict( + "positioner" => fill(eltype(unidi_samps["positioner"])(pmin), pad_nsamps), + "laser" => fill(false, pad_nsamps), + "camera" => fill(false, pad_nsamps), + "stimulus" => fill(false, pad_nsamps) +) + +# We handle padding the traces depending on whether the stimulus signal needs to occur before or after recording starts for a trial +if stim_delay > 0s + append!(pos, "pad_stack_pos", pad_samps["positioner"]) + append!(las1, "pad_stack_las1", pad_samps["laser"]) + append!(cam1, "pad_stack_cam1", pad_samps["camera"]); + prepend!(stim1, "pad_stack_stim1", pad_samps["stimulus"]) #The injection trigger happens after the start of recording +elseif stim_delay < 0s + prepend!(pos, "pad_stack_pos", pad_samps["positioner"]) + prepend!(las1, "pad_stack_las1", pad_samps["laser"]) + prepend!(cam1, "pad_stack_cam1", pad_samps["camera"]) + append!(stim1, "pad_stack_stim1", pad_samps["stimulus"]) #The injection trigger happens before the start of recording +end + +# visualize +# using ImaginePlots, Plots; gr() +# ImaginePlots.plot([pos; las1; cam1; stim1]) + +# save output +fname = "stimulus_control_example.json" +nstacks = sweeps_per_stim * n_trials +nframes_per_stack = Int(sweep_nframes) +write_commands(fname, ocpi2, nstacks, nframes_per_stack, exp_time; isbidi = false) \ No newline at end of file diff --git a/src/ImagineInterface.jl b/src/ImagineInterface.jl index cd952f5..05fffd5 100644 --- a/src/ImagineInterface.jl +++ b/src/ImagineInterface.jl @@ -9,7 +9,7 @@ import ImagineHardware:samprate using AxisArrays const axes = Base.axes -import Base: convert, show, length, size, isempty, ==, append!, pop!, empty!, replace!#, scale +import Base: convert, show, length, size, isempty, ==, append!, prepend!, pop!, empty!, replace!#, scale using Unitful: μm, s, Hz, V diff --git a/src/imaginesignal.jl b/src/imaginesignal.jl index 6a57796..caa5418 100644 --- a/src/imaginesignal.jl +++ b/src/imaginesignal.jl @@ -372,6 +372,49 @@ function append!(com::ImagineSignal{T}, seqname::AbstractString, sequence::T) wh return com end +function prepend!(com::ImagineSignal{T}, seqname::AbstractString) where T<:RLEVector + seqdict = sequence_lookup(com) + if !haskey(seqdict, seqname) + error("The requested sequence name was not found. You most first add the sequence with add_sequence!(com, seqname, sequence), or instead you can add it and append it at the same time with append!(com, seqname, sequence)") + else + #TODO: run safety checks here + #find the length of this sequence and prepend to cumlength vector + seqi = findfirst(x->x==seqname, sequence_names(com)) + pushfirst!(sequence_names(com), seqname) + newseq = sequence_lookup(com)[seqname] + pushfirst!(sequences(com), newseq) + lseq = 0 + if seqi === nothing #we didn't use this sequence yet + lseq = sum(map(count, newseq)) + elseif seqi == 1 + lseq = cumlength(com)[1] + else + lseq = cumlength(com)[seqi] - cumlength(com)[seqi-1] + end + com.cumlength .+= lseq + pushfirst!(com.cumlength, lseq) + end + return com +end + +prepend!(com::ImagineSignal{T}, seqlist::V) where {T<:RLEVector, V<:AbstractVector{String}} = foreach(x->prepend!(com,x), seqlist) + +function prepend!(com::ImagineSignal{T}, seqname::AbstractString, sequence::AbstractVector{TS}) where {T<:RLEVector,TS} + cseq = compress(sequence, mapper(com)) + return prepend!(com, seqname, cseq) +end + +function prepend!(com::ImagineSignal{T}, seqname::AbstractString, sequence::T) where T<:RLEVector + #TODO: run safety checks here + @assert full_length(sequence) >= 1 + add_sequence!(com, seqname, sequence) + pushfirst!(sequences(com), sequence) + pushfirst!(sequence_names(com), seqname) + com.cumlength .+= full_length(sequence) + pushfirst!(com.cumlength, full_length(sequence)) + return com +end + #Repeat the entire sequence currently described by com nreps times #(Equivalent to calling append!(com, seqname) nreps times when seqname is the only sequence in com) function replicate!(com::ImagineSignal{T}, nreps::Int) where T<:RLEVector diff --git a/test/low_level.jl b/test/low_level.jl index 6fef6e4..7dd1ad0 100644 --- a/test/low_level.jl +++ b/test/low_level.jl @@ -293,6 +293,20 @@ append!(pos, "ramp_up2") #append existing, not the first sequence @test length(pos) == lpos + 10 @test all(get_samples(pos, "ramp_up2") .== get_samples(pos)[lpos+6:end]) +#prepend! +empty!(pos; clear_library=true) +newdat = Unitful.μm * [0.0:0.8:800.0...] +append!(pos, "ramp_up", newdat) +stationarydat = Unitful.μm * zeros(50) +prepend!(pos, "stationary", stationarydat) +@test get_samples(pos,1,length(stationarydat)) == stationarydat +prepend!(pos, "stationary") +@test get_samples(pos,1,2*length(stationarydat)) == [stationarydat..., stationarydat...] +dat = get_samples(pos,2*length(stationarydat)+1,length(pos)) +@test dat[1] == mapper(pos).worldmin +@test dat[end] == mapper(pos).worldmax +all([dat[i] > dat[i-1] for i in 2:length(dat)]) + #Test metadata retrieval functions cs2 = chip_size("ocpi-2") cs1 = chip_size("ocpi-1")