Setting up graphs by writing code patterned after stdrun.hoc

[ted]

MOswald wrote:

I use the following command to plot voltage in soma and a dendritic compartment but it the y-axis is set at -1 to 1.

newPlotI()
{
graphItem.addvar("L1CR[0].soma.v(0.5)",2,1)
graphItem.addvar("L1CR[0].dend[65].v(0.5)",3,1)
}

However, if I use the plot command below for the soma only, then the y-axis is scaled automatically.
plotVoltage("Cell.soma.v( 0.5 ),2)

Is there a way to set the y-axis scaling with the first newPlotI procedure?

The bad news is no, but the good news is that it's easy to make your own procedure
that will let you specify whatever y axis scaling you like.

An aside: I'm not aware of any procedure called plotVoltage(). Is that something you
found somewhere in ModelDB or on the WWW? It doesn't seem to be built into NEURON.

Back to the main thread of this discussion--

newPlotI() is the very procedure that is executed when you use the NEURON Main Menu
toolbar to spawn a new "Current axis" graph. This procedure, and its close cousins
newPlotV() and newPlotS(), are defined in stdrun.hoc, which you will find in
nrn/share/nrn/lib/hoc under UNIX / Linux / OS X, or c:\nrnxx\lib\hoc under MSWin (where
the xx is the version number). In fact, here is the code in stdrun.hoc that generates the
"dropdown menu" that appears when you click on the toolbar's Graph button:

Code: Select all

proc graphmenu() {
    xmenu("Graph", 1)
        xbutton("Voltage axis", "newPlotV()")
        xbutton("Current axis", "newPlotI()")
        xbutton("State axis", "newPlotS()")
        xbutton("Shape plot", "newshapeplot()")
        xbutton("Vector movie", "newvectorplot()")
        xbutton("Phase Plane", "newphaseplane()")
        xbutton("Grapher", "load_file(\"grapher.hoc\") makegrapher(1)")
    xmenu()
}

You may want to read about xmenu() and xbuttion() in the Programmer's Reference.

Scrolling a bit farther into stdrun.hoc, we discover the definitions of newPlotV(),
newPlotI(), and newPlotS():

Code: Select all

proc newPlotV() {
        newPlot(0,tstop,-80,40)
        graphItem.save_name("graphList[0].")
        graphList[0].append(graphItem)
        graphItem.addexpr("v(.5)")
}
proc newPlotI() {
        newPlot(0,tstop,-1,1)
        graphItem.save_name("graphList[1].")
        graphList[1].append(graphItem)
}
proc newPlotS() {
        newPlot(0,tstop,0,1)
        graphItem.save_name("graphList[2].")
        graphList[2].append(graphItem)
}

which leads one to guess that the Graphs spawned by any of these three procs will
have an x axis that runs from 0 to tstop, but the y axis range will be
-80 to 40 mV for a "Voltage axis" plot,
-1 to 1 for a "Current axis" plot, and
0 to 1 for a "State axis" plot.
And of course this guess is correct.

So suppose you want a Voltage axis graph, but you want the y axis to span a range
other than -80 to +40? Easy. Write your own procedure that accepts the min and max
y values as arguments. A first stab at this would be

Code: Select all

proc myPlotV() {
        newPlot(0,tstop,$1,$2)
        graphItem.save_name("graphList[0].")
        graphList[0].append(graphItem)
//        graphItem.addexpr("v(.5)")
}

So
myPlotV(-10, 10)
should bring up a graph whose y axis runs from -10 to +10. Try it. Notice that I commented
out graphItem.addexpr("v(.5)"), so this new graph won't come up with membrane
potential at the middle of the default section as its first variable to plot. You'll have to use
graphItem.addvar() to specify what you want plotted (that's probably what you want to
do anyway).

Also notice that I haven't bulletproofed this code at all--what will happen if the first
argument is larger than the second? Nonsense results can be prevented by inserting a
conditional statement just before the newPlot() call, e.g.

Code: Select all

proc myPlotV() {
        if ($1>=$2) {
          print "myPlotV error--first argument must be smaller than second argument"
          stop
        }
        newPlot(0,tstop,$1,$2)
        graphItem.save_name("graphList[0].")
        graphList[0].append(graphItem)
//        graphItem.addexpr("v(.5)")
}

Why are there different graphLists?

Graphs spawned by newPlotV(), newPlotI(), and newPlotS() are appended to one of
three graphLists, so that they will be automatically updated during simulation execution.
As noted at the end of chapter 7 in The NEURON Book, when fixed time step integration
is used, graphList[0], [1], and [2] differ in the x coordinates used for plotted points. For
graphList[0] each variable is plotted vs. the current value of t, but for graphList[1] the x
coordinates are t - 0.5*dt, and for graphList[2] the x coords are t + 0.5*dt.

The reason for this is that, if you use NEURON's Crank-Nicholson integrator (by
assigning secondorder = 2 -- see the Programmer's Reference documentation of
secondorder), voltages will be second order correct at integral multiples of dt, but
currents and states will be second order correct half a time step earlier or later,
respectively. By plotting voltages in Voltage axis graphs, currents in Current axis graphs,
and states in State axis graphs, the plotted points preserve second order accuracy in
time.

"But doesn't this mess up the accuracy of current and state graphs when the default
backward Euler method is used?" No, because backward Euler is only first order
accurate in time, i.e. local error of the solution is proportional to dt, the size of the time
step. This means that there isn't a specific time in the interval from (t - 0.5*dt) to (t + 0.5*dt)
at which voltages, currents, and states are first order correct--they're equally
(in)accurate over that entire interval.

When adaptive integration ("CVODE") is used, all variables are computed to the same
degree of accuracy at identical times, and all "Voltage axis," "Current axis," and "State
axis" graphs plot their variables at identical values of t.

[MOswald]

Thank you for this detailed answer to my plotting question.

I can see that I should be using the newPlotV() command rather than that for current, and how to specify the axis limits if I want to change them from the default axis settings.

The plotVoltage() command below was provided to me by another Neuron computational group as part of a hoc file that specifies some common procedures used in a typical simulation. At least I understand this better now :)
Cheers, Manfred

/*****************************************************************

This function is used for voltage plot

Required Inputs : $s1. what to plot

$2. color

*****************************************************************/
proc plotVoltage(){ localobj grph
grph = new Graph(0)
//grph.addexpr($s1, 2, 1, 0.541214, 1.04377, 2)
if(numarg() == 1){
//grph.addexpr($s1, 1, 1 )
grph.addexpr("",1,1)
}else{
//grph.addexpr($s1, $2, 1 )
grph.addexpr("",$s1,$2,1)
}
graphList[0].append(grph)
grph.view(0, -80, tstop, 120,300, 75, 800, 325)
grph.label(0.4,0.9,$s1)
}

So looking up the SYNTAX for .view function I figure out that this does also set coordinates for the graph axis as well as setting the window coordinates relative to the display field overall:

.view(mleft, mbottom, mwidth, mheight, wleft,
wtop, wwidth, wheight)
.view(2)
DESCRIPTION
Map a view of the Shape scene. m stands for model coordinates within the window, w stands for screen coordinates for placement and size of the window. The placement of the window with respect to the screen is intended to be precise and is with respect to pixel coordinates where 0,0 is the top left corner of the screen.

The single argument form maps a view in which the aspect ratio between x and y axes is always 1. eg like a shape window.

[MOswald]

PS. Could I modify the plotVoltage() procedure to plot two or more variables into the same graph window?

I have tried addvar and addexpr commands in various ways but without success so far.

Ta, Manfred

[ted]

Only at the cost of increased complexity that would impede readability and maintainability of
code. Your best bet is to add variables with individual statements outside of a proc, taking
advantage of the fact that the most recently created graph is the last item in the graphList.

[ted]

I will amend my previous remark slightly. Here's a first stab at something
that would not be too onerous to implement and use. Instead of passing
a variable number of arguments, pass your proc a List whose contents are
the seeds of each command that you want addvar to execute. Then iterate
through the List, using sprint to create each command string, and executing
that command string, before moving on to the next command seed in the List.

Code: Select all

objref argList, argStr
argList = new List()
tmpStr = new String() // an object wrapper for strdef--see begintemplate String 
  // in stdlib.hoc, which is automatically loaded if you run nrngui, 
  // or execute load_file("nrngui.hoc") or load_file("noload.hoc")
// before calling your proc, fill the List with String objects, each of which contains 
// a strdef that is part of the argument that you eventually want to use in 
// an addexpr call
// An example that uses bogus variable names:
tmpStr.s = "\"foo\", 1" // name in quotes and color
argList.append(tmpStr)
tmpStr.s = "\"bah\", 2"
argList.append(tmpStr)
tmpStr.s = "\"woo\", 3"
argList.append(tmpStr)

// expects a List of String objrefs as its only arg
strdef tmpstr
proc myplotv() { local i  localobj g
  g = new Graph(0)
  for i = 0, $o1.count()-1 {
    sprint(tmpstr, "g.addexpr(%s, 1)", $o1.object(i).s) ) // this example uses brush 1 for all
    print(tmpstr) // check the command strings before you comment this out
      // and uncomment the next statement
//    execute(tmpstr)
  }
  graphList[0].append(g)
  // you may wish to replace some of these magic numbers
  g.view(0, -80, tstop, 120,300, 75, 800, 325)
  // with symbolic "constants" whose values are defined elsewhere
  // or even passed as arguments to this proc, e.g.
//  g.view(0, VMIN, tstop, VMAX-VMIN, 300, 75, 800, 325)
  // ditto for graph position and dimensions
}