/*
* $Id: stim.hoc,v 1.5 2009/02/24 00:55:27 ted Exp ted $
* Set up stimulation waveform
* 2018/05/20 Modified by Aman Aberra
*/

// default values
DEL = 1  // ms - delay to first phase
DUR = 1 // ms - duration of first phase
AMP = -10  // µA (stim_mode = 1) or V/m (stim_mode = 2) - amplitude of first phase
objref stim_amp, stim_time
stim_amp = new Vector()
stim_time = new Vector()

proc stim_waveform() {
  // this uses interpolated play
  // index    0  1    2    3        4        5          
  // stim vec 0, 0,   1,   1,       0        0          
  // time vec 0, DEL, DEL, DEL+DUR, DEL+DUR, DEL+DUR+1
  //  really  0, $1,  $1,  $1+$2,   $1+$2,   $1+$2+1
  // first the stim vector
  stim_amp.resize(6)
  stim_amp.fill(0)
  stim_amp.x[2]=1
  stim_amp.x[3]=1  
  stim_amp.mul($3)
  // now the time vector
  stim_time.resize(6)
  stim_time.x[1]=$1
  stim_time.x[2]=$1
  stim_time.x[3]=$1+$2
  stim_time.x[4]=$1+$2
  stim_time.x[5]=$1+$2+1
}


ATTACHED__ = 0

proc attach_stim() {
// since stim_xtra is GLOBAL, we only need to specify Vector.play()
// for one instance of xtra, i.e. at just one internal node
// of only one section that contains xtra
  forall {  // check each section to find one that has xtra
    if (ATTACHED__ == 0) {  // don't bother if stim is already attached to something
      if (ismembrane("xtra")) {
        stim_amp.play(&stim_xtra, stim_time, 1) // "interpolated" play
        ATTACHED__ = 1
      }
    }
  }
}

objref g1
proc plot_waveform() {
	g1 = new Graph(0)
	g1.size(0,stim_time.size(),-1, $1)
	stim_amp.plot(g1,stim_time)
  if ($1 > 0){ // anodic pulse
    g1.view(0,-$1/2,stim_time.x[stim_time.size()-1],$1*2,1081, 547, 300.48, 200.32) 
  }	else { // cathodic pulse
    g1.view(0,$1*3/2,stim_time.x[stim_time.size()-1],-$1*2,1081, 547, 300.48, 200.32) 
  }
}

proc setstim() {
  del = $1
  dur = $2
  amp = $3
  stim_waveform(del, dur, amp)
  attach_stim()
  plot_waveform(amp)
  printf("Generated waveform with del = %g ms, dur = %g ms, amp = %g uA or V/m\n",del,dur,amp)		
}

xpanel("Temporal parameters for extracellular stimulation", 0)
  xvalue("Delay (ms)", "DEL", 1, "setstim(DEL,DUR,AMP)", 0, 1)
  xvalue("Duration (ms)", "DUR", 1, "setstim(DEL,DUR,AMP)", 0, 1)
  xvalue("Amplitude (uA or V/m)", "AMP", 1, "setstim(DEL,DUR,AMP)", 0, 1)
xpanel(535,652)