-
Notifications
You must be signed in to change notification settings - Fork 2
/
Copy pathdecvec.hoc
451 lines (406 loc) · 11.8 KB
/
decvec.hoc
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
// $Id: decvec.hoc,v 1.63 2001/06/14 21:13:40 billl Exp $
proc decvec() {}
//* Declarations
objref ind, vec, vec0, vec1, tmpvec, vrtmp, veclist
objref tmpobj, XO, YO, rdm
load_file("declist.hoc") // declare list iterators
print "Loading decvec"
{ MSONUM=100 MSOSIZ=100 msomax=0 msoptr=0 objref mso[MSONUM] }
double x[4],y[4]
xx=0 // declare a scalar
ind = new Vector(100)
vec = new Vector(100)
vec0 = new Vector()
vec1 = new Vector()
vrtmp = new Vector()
veclist = new List()
rdm = new Random()
if (! xwindows) {
objref graphItem
strdef temp_string_, temp_string2_
}
strdef xtmp
if (wopen("xtmp")) xtmp = "xtmp" else xtmp="/tmp/xtmp" // scratch file to save system output to
//* stuff that doesn't belong here
objref dir
dir = new List()
//** dired(list,file) - put together list of files matching 'file', calls 'ls -1 file'
// dired(list,file,1) file name to read for list of files
proc dired () {
if (numarg()==0) { print "dired(list,filename) adds the filename to list (use wildcards)"
return }
if (numarg()==3) {
tmpfile.ropen($s2)
} else {
sprint(temp_string_,"ls -1R %s > %s",$s2,xtmp) // list in order of creation time
system(temp_string_)
tmpfile.ropen(xtmp)
}
while (tmpfile.scanstr(temp_string_) != -1) {
tmpobj=new String()
tmpobj.s=temp_string_
$o1.append(tmpobj)
tmpfile.gets(temp_string_) // get rid of the rest of the line
}
}
//** lbrw(list,name,action) is used to put up a browser
// note action given without '()'
proc lbrw () {
$o1.browser($s2,"s")
sprint($s2,"%s()",$s2)
$o1.accept_action($s2)
}
//* vector iterator vtr
// usage 'for vtr(&x, vec) { print x }'
iterator vtr() { local i
if (numarg()==3) {$&3=0} else {i1 = 0}
for i = 0, $o2.size() - 1 {
$&1 = $o2.x[i]
iterator_statement
if (numarg()==3) { $&3+=1 } else { i1+=1 }
}
}
//* vector iterator vtr2, treat two vectors as pairs
// usage 'for vtr2(&x, &y, vec1, vec2) { print x,y }'
iterator vtr2() { local i
if ($o3.size != $o4.size) { print "vtr2 ERROR: sizes differ." return }
if (numarg()==5) {$&5=0} else {i1 = 0}
for i = 0, $o3.size() - 1 {
$&1 = $o3.x[i]
$&2 = $o4.x[i]
iterator_statement
if (numarg()==5) { $&5+=1 } else { i1+=1 }
}
}
//* iterator lvtr, step through a list and a vector together
// usage 'for lvtr(XO, &x, list, vec) { print XO,x }'
iterator lvtr() { local i
if ($o3.count < $o4.size) { print "lvtr ERROR: vecsize>listsize." return }
if ($o3.count != $o4.size) { print "lvtr WARNING: sizes differ." }
if (numarg()==5) {$&5=0} else {i1 = 0}
for i = 0, $o4.size()-1 {
$o1 = $o3.object(i)
$&2 = $o4.x[i]
iterator_statement
if (numarg()==5) { $&5+=1 } else { i1+=1 }
}
}
//* other iterators: case, scase, ocase
iterator case() { local i
i1 = 0
for i = 2, numarg() {
$&1 = $i
iterator_statement
i1+=1
}
}
iterator scase() { local i
i1 = 0
for i = 1, numarg() {
temp_string_ = $si
iterator_statement
i1+=1
}
}
// eg for scase2("a","b","c","d","e","f") print temp_string_,temp_string2_
iterator scase2() { local i
i1 = 0
if (numarg()%2==1) {print "ERROR: scase2 needs even number of args" return }
for i = 1, numarg() {
tmpobj=new String2()
tmpobj.s=$si i+=1 tmpobj.t=$si
iterator_statement
i1+=1
}
}
iterator ocase() { local i
i1 = 0
for i = 1, numarg() {
XO = $oi
iterator_statement
}
i1+=1
}
//* nind(targ,data,ind) fill the target vector with data NOT index by ind (opposite of v.index)
proc nind () {
if (! eqobj($o1,$o2)) $o1.copy($o2)
$o1.indset($o3,-1e20)
$o1.where($o1,">",-1e20)
}
//* vlk(vec)
// vlk(vec,max)
// vlk(vec,min,max)
// prints out a segment of a vector
vlk_width=20
proc vlk () { local i,j,min,max,dual,wdh,nonl
j=dual=0 nl=1 wdh=vlk_width
if (numarg()==1) { min=0 max=$o1.size-1 }
if (numarg()==2) if ($2==0) { nl=min=0 max=$o1.size-1 // vlk(vec,0) flag to suppress new lines
} else if ($2>0) { min=0 max=$2-1 } else { min=$o1.size+$2 max=$o1.size-1 }
if (numarg()==3) if ($3>-1) { min=$2 max=$3 } else { min=0 max=$o1.size-1 dual=1 }
if (numarg()==4) { min=$3 max=$4 dual=1 }
if (min<0) min=0
if (max>$o1.size-1) max=$o1.size-1
if (dual) if (max>$o2.size-1) max=$o2.size-1
for i=min,max {
if (dual) printf("%g:%g ",$o1.x[i],$o2.x[i]) else printf("%g ",$o1.x[i])
if ((j=j+1)%vlk_width==0 && nl) { print "" }
}
if (nl) print ""
}
//* vprf() prints 1,2 or 3 vectors in parallel to output file
proc vprf () { local x2
if (! tmpfile.isopen()) {
print "Writing to temp file 'temp'"
tmpfile.wopen("temp")
}
if (numarg()==1) {
for vtr(&x,$o1) { tmpfile.printf("%g\n",x) }
} else if (numarg()==2) {
for vtr2(&x,&y,$o1,$o2) { tmpfile.printf("%g %g\n",x,y) }
} else if (numarg()==3) {
for vtr2(&x,&y,$o1,$o2,&ii) { x2=$o3.x[ii] tmpfile.printf("%g %g %g\n",x,y,x2) }
}
tmpfile.close
}
//* vpr() prints 1,2 or 3 vectors in parallel to STDOUT
proc vpr () { local x2
if (numarg()==1) {
for vtr(&x,$o1) { printf("%g ",x) }
} else if (numarg()==2) {
for vtr2(&x,&y,$o1,$o2) { printf("%g:%g ",x,y) }
} else if (numarg()==3) {
for vtr2(&x,&y,$o1,$o2,&ii) { x2=$o3.x[ii] printf("%g:%g:%g ",x,y,x2) }
}
print ""
}
//* readvec(vec) read from line
proc readvec () {
$o1.resize(0)
while (read(xx)) $o1.append(xx)
vlk($o1)
}
//* popvec(), savenums, readnums, vecsprint, savevec
proc pushvec () { local i // same as .append, retained for compatability
for i=2, numarg() $o1.append($i)
}
proc revec () { local i // clear vector then append
$o1.resize(0)
for i=2, numarg() $o1.append($i)
}
// savenums(x[,y,...]) save numbers to tmpfile via a vector
proc savenums () { local vv,i
vv=allocvecs(1)
for i=1, numarg() mso[vv].append($i)
mso[vv].vwrite(tmpfile)
dealloc(vv)
}
// readnums(&x[,&y...]) recover nums from tmpfile via a vector
proc readnums () { local vv,i
vv=allocvecs(1)
mso[vv].vread(tmpfile)
if (numarg()!=mso[vv].size) { print "Error: size in readnums" return }
for i=1,numarg() $&i = mso[vv].x[i-1]
dealloc(vv)
}
//* remove last entry
func popvec () { local sz, ret
sz = $o1.size-1
ret = $o1.x[sz]
$o1.resize[sz]
return ret
}
//* chkvec (look at last entry)
func chkvec () { if ($o1.size>0) return $o1.x[$o1.size-1] else return -1e10 }
// vecsprint(strdef,vec)
proc vecsprint () { local ii
if ($o2.size>100) { return }
for ii=0,$o2.size-1 { sprint($s1,"%s %g ",$s1,$o2.x[ii]) }
}
// savevec([list,]vec1[,vec2,...]) add vector onto veclist or other list if given as 1st arg
proc savevec () { local i,flag,beg
if (isobj($o1,"List")) beg=2 else beg=1
for i=beg, numarg() {
tmpvec = new Vector($oi.size)
tmpvec.copy($oi)
if (beg==2) $o1.append(tmpvec) else veclist.append(tmpvec)
tmpvec = nil
}
}
// tvecl() -- transpose veclist
proc tvecl () { local cnt,sz,err,ii,p
err = 0
cnt = veclist.count
sz = veclist.object(0).size
for ltr(XO,veclist) if (XO.size!=sz) err=i1
if (err) { print "Wrong size vector is #",i1 return }
p = allocvecs(1,cnt) mso[p].resize(cnt)
for ii=0,sz-1 {
for jj=0,cnt-1 {
XO=veclist.object(jj)
mso[p].x[jj] = XO.x[ii]
}
savevec(mso[p])
}
for (jj=cnt-1;jj>=0;jj-=1) veclist.remove(jj)
}
//* vinsect(v1,v2,v3) -- v1 gets intersection (common members) of v2,v3
// replaced by v.insct() in vecst.mod
proc vinsect () {
$o1.resize(0)
for vtr(&x,$o2) for vtr(&y,$o3) if (x==y) $o1.append(x)
}
//* vecsplit(vec,vec1,vec2[,vec3,...])
// splits vec into other vecs given
proc vecsplit () { local num,ii,i
num = numarg()-1 // how many
for i=2,numarg() $oi.resize(0)
for (ii=0;ii<$o1.size;ii+=num) {
for i=2,numarg() if (ii+i-2<$o1.size) $oi.append($o1.x[ii+i-2])
}
}
//* vecsort(vec,vec1,vec2[,vec3,...])
// sorts n vecs including first vec by first one
proc vecsort () { local i,inv,scr,narg
narg=numarg()
if (narg<2 || narg>10) {print "Wrong #args in decvec.hoc:vecsort" return}
scr=inv=allocvecs(2) scr+=1
$o1.sortindex(mso[inv])
mso[scr].resize(mso[inv].size)
sprint(temp_string_,"%s.fewind(%s,%s,%s",mso[scr],mso[inv],$o1,$o2)
for i=3,narg sprint(temp_string_,"%s,%s",temp_string_,$oi)
sprint(temp_string_,"%s)",temp_string_)
execute(temp_string_)
dealloc(inv)
}
//* vdelind() -- delete a single index
proc vdelind () { local i,iin
iin = $2
if (iin<0) iin=$o1.size+iin
if (iin>$o1.size-1 || iin<0) {
printf("vdelind Error: index %d doesn't exist.\n",iin) return }
if (iin<$o1.size-1) $o1.copy($o1,iin,iin+1,$o1.size-1)
$o1.resize($o1.size-1)
}
//* mkveclist(num[,sz]) recreate veclist to have NUM vecs each of size SZ (or MSOSIZ)
proc mkveclist () { local ii,num,sz,diff
num=$1
diff = num-veclist.count
if (numarg()==2) { sz=$2 } else { sz = MSOSIZ }
if (diff>0) {
for ii=0,diff-1 {
tmpvec = new Vector(sz)
veclist.append(tmpvec)
}
} else if (diff<0) {
for (ii=veclist.count-1;ii>=num;ii=ii-1) { veclist.remove(ii) }
}
for ltr(XO,veclist) { XO.resize(sz) }
}
//* allocvecs
// create temp set of vectors on mso
// returns starting point on mso
// eg p = allocvecs(3)
// access these vectors by mso[p+0] ... [p+2]
func allocvecs () { local ii, llen, sz, newv
if (numarg()==0) { print "p=allocvecs(#), access with mso[p], mso[p+1]..." return 0}
newv = $1
if (numarg()==2) { sz=$2 } else { sz=MSOSIZ }
llen = msoptr
for ii=msomax,msoptr+newv-1 { // may need new vectors
if (ii>=MSONUM) { print "alloc ERROR: MSONUM exceeded." return 0 }
mso[ii] = new Vector(sz)
}
for ii=0,newv-1 {
mso[msoptr].resize(0)
msoptr = msoptr+1
}
if (msomax<msoptr) msomax = msoptr
return llen
}
//** dealloc(start)
// remove temp set of vectors from mso
proc dealloc () { local ii,min
if (numarg()==0) { min = 0 } else { min = $1 }
msomax = msoptr
msoptr = min
}
//* vecconcat(vec1,vec2,...)
// destructive: concatenates all vecs onto vec1
proc vecconcat () { local i
if (numarg()<2) { print "vecconcat(v1,v2,...) puts all into v1" return }
for i=2,numarg() {
$o1.copy($oi,$o1.size)
}
}
//** vecelim(v1,v2) eliminates items in v1 given by index vec v2
proc vecelim () {
for vtr(&x,$o2) { $o1.x[x]= -1e20 }
$o1.where($o1,"!=",-1e20)
}
//** veceq() like vec.eq but don't have to be same size and shows discrepency
func veceq () { local sz1,sz2,eq,beg,ii,jj,kk
sz1=$o1.size sz2=$o2.size
if (numarg()==3) beg=$3 else beg=0
if (sz1!=sz2) printf("%s %d; %s %d\n",$o1,sz1,$o2,sz2)
ii=0 jj=beg
while (ii<sz1 && jj<sz2) {
if ($o1.x[ii]!=$o2.x[jj]) {
eq=0
printf("Differ at %d %d\n",ii,jj)
for kk=-10,10 if ((ii+kk)>=0 && (ii+kk)<sz1 && (jj+kk)>=0 && (jj+kk)<sz2) {
printf("(%d)%g:(%d)%g ",(ii+kk),$o1.x[ii+kk],(jj+kk),$o2.x[jj+kk]) }
print ""
break
} else eq=1
ii+=1 jj=ii+beg
}
return eq
}
//* isstring() determine if object $o1 is of type string, if so return the string in $s2
func isstring () {
sprint($s2,"%s",$o1)
if (sfunc.substr($s2,"String")==0) {
sprint($s2,"%s",$o1.s)
return 1
} else {
sprint($s2,"%s",$o1)
return 0
}
}
//** isassigned() checks whether an object is Null
func isassigned () {
sprint(temp_string_,"%s",$o1)
if (sfunc.substr(temp_string_,"NULLobject")==0) {
return 0
} else {
return 1
}
}
//** eqobj(o1,o2) checks whether 2 objects are the same
func eqobj () { return object_id($o1) == object_id($o2) }
//** isobj(o1,s2) checks whether object $o1 is of type $s2
func isobj () {
sprint(temp_string_,"%s",$o1)
if (sfunc.substr(temp_string_,$s2)==0) {
return 1
} else {
return 0
}
}
// destructive of $s1
func str2num () { local ii
sscanf($s1,"%d",&x)
return x
}
// like perl chop -- removes the last character
proc chop () { sfunc.left($s1,sfunc.len($s1)-1) }
// newlst() puts a newline in the middle of a string
proc newlst () { local l
if (numarg()>1) l=$2 else l=int(sfunc.len($s1)/2)
temp_string_=$s1
temp_string2_=$s1
sfunc.left(temp_string_,l)
sfunc.right(temp_string2_,l)
sprint($s1,"%s\n%s",temp_string_,temp_string2_)
}