Page 1 of 1

how to use savestate() with user-defined point processes

Posted: Fri Oct 17, 2008 12:30 pm
by sabsan
Dear all,
there is something that I've recently noted and I cannot explain. I wrote the following .mod file about a point process that may generate a train of square pulses:

Code: Select all

 :Current Clamp

 NEURON {
 		 POINT_PROCESS trainIClamp
 		 RANGE del, PW, train, amp, freq, i, conv, pulsecount, onoff
 		 ELECTRODE_CURRENT i
 }

 UNITS { (na) = (nanoamp) }

 PARAMETER{
 		del (ms)
 		PW (ms)
 		train (ms)
 		amp (na)
 		freq (1/s)
 		conv = 1000 (ms/s)
 		pulsecount (s/s)
 		onoff (s/s)
 }

 ASSIGNED {
 		i (na)		
 }

 INITIAL  { LOCAL j,k
 			pulsecount = 0
 			onoff = 0
            k =  (train/conv)/freq
 			i = 0
 			FROM j = 0 TO k  {
 				at_time (del + (j*(conv/freq)))
		 		at_time (del + PW + (j*(conv/freq)))
 		  	}
 		  	at_time (del + train)
 }

 BREAKPOINT {
		 		if (t < del + train && t > del) {
		 				if (t > del + (pulsecount*(conv/freq)) && t < del + (pulsecount*(conv/freq)) + PW)  {
		 						i = amp
		 						onoff = 1
		 				} else {
		 						if (onoff == 0) {
		 							i = 0
		 						} else {
		 							i = 0
		 							pulsecount = pulsecount + 1
		 							onoff = 0
		 						}
		 				}
		 		} else {
		 				i = 0
		 				pulsecount = 0
		 				onoff = 0
		 		}
 		
 }
User can customize amplitude, frequency, initial delay, pulsewidth and duration of the train. I used this point process in a dummy example (I inserted it in a single compartment soma). I run the example for a while (500 ms) and saved the state in a SaveState object and stored it in a file (if interested, please, execute file "test1.hoc" at the end of this post). Then, I run a second file ("test2.hoc" at the end of this post), restoring the previously saved state, and noted that the pulses of my train are not delivered anymore, even if the features of the train in both test1 and test2 file are the same. I guess that the .restore() command clears up the events queued by the at_time() commands of my point process. Is it correct? Do you have any suggestion to avoid it? Thank you very much. Bye.

test1.hoc:

Code: Select all

// test no.1: save the current state

load_file("nrngui.hoc")

create soma

soma {
	nseg = 1
    diam = 18.8
    L = 18.8
    Ra = 123.0
    
    insert hh
    ena = 71.5
    ek = -89.1
    gnabar_hh=0.25
    gl_hh = .0001666
    el_hh = -60.0
}

objref inIClmp, f, ss
    
soma {
	inIClmp = new trainIClamp(0.5)
    inIClmp.del = 200
    inIClmp.PW = 2
    inIClmp.train = 1e9
    inIClmp.freq = 50
    inIClmp.amp = 1
}


f = new File("state.01")
ss = new SaveState()

proc savestate() {
	ss.save()
	ss.fwrite(f)
}

load_file("sthD1.ses")
finitialize(-80)
fcurrent()

tstop=500
run()
savestate()
f.close()
print "The end"
test2.hoc:

Code: Select all

// test no.2: restore the saved state

load_file("nrngui.hoc")

create soma

soma {
	nseg = 1
    diam = 18.8
    L = 18.8
    Ra = 123.0
    
    insert hh
    ena = 71.5
    ek = -89.1
    gnabar_hh=0.25
    gl_hh = .0001666
    el_hh = -60.0
}

objref inIClmp, f, ss
    
soma {
	inIClmp = new trainIClamp(0.5)
    inIClmp.del = 200
    inIClmp.PW = 2
    inIClmp.train = 1e9
    inIClmp.freq = 50
    inIClmp.amp = 1
}


f = new File("state.01")
ss = new SaveState()

proc init() {
	ss.fread(f,0)
	finitialize(-80)
	ss.restore()
  	
	if (cvode.active()) {
		cvode.re_init()
	} else {
		fcurrent()
	}
	frecord_init()
}

load_file("sthD1.ses")

init()
tstop=800
run()
f.close()
print "The end"
A .ses file to plot the train of pulses (it is named "sthD1.ses"):

Code: Select all

{load_file("nrngui.hoc")}
objectvar save_window_, rvp_
objectvar scene_vector_[4]
objectvar ocbox_, ocbox_list_, scene_, scene_list_
{ocbox_list_ = new List()  scene_list_ = new List()}
{pwman_place(0,0,0)}
{
save_window_ = new Graph(0)
save_window_.size(-10,800,0,1)
scene_vector_[2] = save_window_
{save_window_.view(-10, 0, 810, 1, 390, 18, 300.6, 200.8)}
graphList[1].append(save_window_)
save_window_.save_name("graphList[1].")
save_window_.addvar("inIClmp.i", 1, 1, 0.8, 0.9, 2)
}
{
xpanel("RunControl", 0)
v_init = -65
xvalue("Init","v_init", 1,"stdinit()", 1, 1 )
xbutton("Init & Run","run()")
xbutton("Stop","stoprun=1")
runStopAt = 5
xvalue("Continue til","runStopAt", 1,"{continuerun(runStopAt) stoprun=1}", 1, 1 )
runStopIn = 1
xvalue("Continue for","runStopIn", 1,"{continuerun(t + runStopIn) stoprun=1}", 1, 1 )
xbutton("Single Step","steprun()")
t = 800
xvalue("t","t", 2 )
tstop = 800
xvalue("Tstop","tstop", 1,"tstop_changed()", 0, 1 )
dt = 0.025
xvalue("dt","dt", 1,"setdt()", 0, 1 )
steps_per_ms = 40
xvalue("Points plotted/ms","steps_per_ms", 1,"setdt()", 0, 1 )
screen_update_invl = 0.05
xvalue("Scrn update invl","screen_update_invl", 1,"", 0, 1 )
realtime = 0.47
xvalue("Real Time","realtime", 0,"", 0, 1 )
xpanel(823,114)
}
{
save_window_ = new Graph(0)
save_window_.size(0,800,-90,70)
scene_vector_[3] = save_window_
{save_window_.view(0, -90, 800, 160, 54, 356, 300.6, 200.8)}
graphList[0].append(save_window_)
save_window_.save_name("graphList[0].")
save_window_.addvar("soma.v( 0.5 )", 1, 1, 0.8, 0.9, 2)
}
objectvar scene_vector_[1]
{doNotify()}

Re: how to use savestate() with user-defined point processes

Posted: Sat Oct 18, 2008 12:19 pm
by ted
Two comments, neither of which will fix the problem you encountered--

1. The "conv" variable should not be RANGE, unless you really want to allow the possibility that a hoc statement or careless use of a GUI tool will radically alter the pules times and durations. It is simply a scale factor, and should be replaced by the numeric value 1000, which you may enclose in parentheses so that modlunit will not complain. For example,
at_time (del + (j*(conv/freq)))
should be written
at_time (del + (j*((1000)/freq)))

2. at_time() is deprecated. It still works, but it should not be used in new code. The preferred way to do this stuff is with events, which generally results in code that is easier to write and debug. For two mod files that use events to generate trains of current pulses, see
periodic stimulation?
viewtopic.php?f=8&t=137&p=275

That said, I tried saving and restoring states with a model that used IPulse2, and I'm seeing segmentation errors. So either there's a bug in my rather simple hoc code, or a bug in the state save/restore functions with the version of NEURON that I'm running (7.0 (171:27d1e43054cc) 2008-09-26)

Re: how to use savestate() with user-defined point processes

Posted: Sun Oct 19, 2008 7:16 am
by hines
That is a bug in SaveState. Probably happens when a SaveState.read is done
and a net_send event from an INITIAL block is still on the queue.
Also SaveState is not yet functional in the thread (7.0) version.

SaveState alwasy seems to be the last thing updated for any extension.
I'll work on making it work again.

Re: how to use savestate() with user-defined point processes

Posted: Sun Oct 19, 2008 12:36 pm
by ted
Ipulse2 does work properly with NEURON 6.2--no segmentation error.

SaveState.restore() sets up a simulation in which Ipulse2 produces only those current pulses that were not already generated in the run whose states were saved--even if you execute t = 0 right after the restore() statement. Why? Because the custom proc init() calls ss.restore() after it calls finitialize()--

Code: Select all

proc init() {
  finitialize(v_init)
  ss.restore()
//  t = 0 // won't affect events

  if (cvode.active()) {
    cvode.re_init()
  } else {
    fcurrent()
  }
  frecord_init()
}
The SaveState documentation http://www.neuron.yale.edu/neuron/stati ... state.html says
The state also includes all the outstanding events (external and self) . . .

The outstanding event delivery times are absolute. When restored, all outstanding events will be cleared and the restored event times and NetCon info will take their place.
"Outstanding events" means events that have been launched but have not yet been been delivered.

Re: how to use savestate() with user-defined point processes

Posted: Mon Oct 20, 2008 4:17 pm
by hines
SaveState functionality has been restored in the 7.0 (thread) version.
http://www.neuron.yale.edu/hg/neuron/nr ... d5c8983552
I have not
remade the alpha distribution files so presently it is available only from the
mercurial repository via
hg clone http://www.neuron.yale.edu/hg/neuron/nrn