// based on interpxyz.hoc,v 1.2 2005/09/10 23:02:15
/* Computes xyz coords of nodes in a model cell 
   whose topology & geometry are defined by pt3d data.
   Expects sections to already exist, and that the xtrau mechanism has been inserted
 */


// original data, irregularly spaced
objref xx, yy, zz, length
// interpolated data, spaced at regular intervals
objref xint, yint, zint, range

proc grindaway() { local ii, nn, kk, xr
  forall {
//    if (ismembrane("xtra")) {
    if (ismembrane("xtrau")) {
    // get the data for the section
      nn = n3d()
      xx = new Vector(nn)
      yy = new Vector(nn)
      zz = new Vector(nn)
      length = new Vector(nn)

      for ii = 0,nn-1 {
        xx.x[ii] = x3d(ii)
        yy.x[ii] = y3d(ii)
        zz.x[ii] = z3d(ii)
        length.x[ii] = arc3d(ii)
      }

      // to use Vector class's .interpolate() 
      // must first scale the independent variable
      // i.e. normalize length along centroid
      length.div(length.x[nn-1])

      // initialize the destination "independent" vector
      range = new Vector(nseg+2)
      range.indgen(1/nseg)
      range.sub(1/(2*nseg))
      range.x[0]=0
      range.x[nseg+1]=1

      // length contains the normalized distances of the pt3d points 
      // along the centroid of the section.  These are spaced at 
      // irregular intervals.
      // range contains the normalized distances of the nodes along the 
      // centroid of the section.  These are spaced at regular intervals.
      // Ready to interpolate.

      xint = new Vector(nseg+2)
      yint = new Vector(nseg+2)
      zint = new Vector(nseg+2)
      xint.interpolate(range, length, xx)
      yint.interpolate(range, length, yy)
      zint.interpolate(range, length, zz)

//      // for each node, assign the xyz values to x_xtra, y_xtra, z_xtra
      // for each node, assign the xyz values to x_xtrau, y_xtrau, z_xtrau
//      for ii = 0, nseg+1 {
// don't bother computing coords of the 0 and 1 ends
// also avoid writing coords of the 1 end into the last internal node's coords
      for ii = 1, nseg {
        xr = range.x[ii]
/*
        x_xtra(xr) = xint.x[ii]
        y_xtra(xr) = yint.x[ii]
        z_xtra(xr) = zint.x[ii]
*/
        x_xtrau(xr) = xint.x[ii]
        y_xtrau(xr) = yint.x[ii]
        z_xtrau(xr) = zint.x[ii]
      }
    }
  }
}