""

svn path=/trunk/; revision=468

25 files changed, 1015 insertions, 839 deletions

diff --git a/externals/grill/vasp/changes.txt b/externals/grill/vasp/changes.txt
index 2baf24ac..0f6e4736 100644
--- a/externals/grill/vasp/changes.txt
+++ b/externals/grill/vasp/changes.txt
@@ -1,119 +1,123 @@
-Version history:

-

-0.1.2:

-- FIX: bug in vasp.frames* ... wrong argument

-- FIX: another bug in resize of immediate buffers

-- CHANGE: immediates will be initialized to zero, also upon resize to a larger size (by default, later make this optional)

-- ADD: attributes for all basic parameters: ref,to,arg, update,detach,prior (loglvl,argchk,unit later)

-- ADD: vasp.channel,vasp.vector - attribute "index"

-- ADD: vasp.frames,vasp.frames+,vasp.offset,vasp.offset+,vasp.size,vasp.size+,vasp.imm - attribute "frames"

-- ADD: vasp.frames*,vasp.frames/,vasp.size*,vasp.size/ - attribute "factor"

-- ADD: vasp.part - attribute "parts"

-- ADD: vasp.peaks? - attribute "peaks"

-- ADD: vasp.shift,vasp.xshift - attribute "fill"

-- ADD: vasp.tilt,vasp.xtilt - attributes "fill","inter"

-- ADD: vasp.offset=,vasp.frames= - attributes "dir","slope"

-- CHANGE: update m_method_ virtual functions of same objects

-- FIX: immediate vasp garbage collection is now done in separate thread

-- ADD: vasp.size... objects: attribute "keep" determines whether buffer content should be preserved (default = 1)

-- CHANGE: buffer reference checks upon all operations enforced (by new CVasp (checked vasp) class)

-- FIX: immediate vasps can not have frames > size

-- FIX: immediate garbage collection was not triggered in threaded build

-- FIX: changed forgotten standard ASSERT to FLEXT_ASSERT (thanks to Frank Barknecht)

-

-0.1.1:

-- FIX: bug in resize of immediate buffers

-- CHANGE: check (at actual operation!) whether vector exceeds buffer size (GetRVecs, GetCVecs -> corrlen)

-

-0.1.0:

-- NEW: vasp.a<, vasp.a>, vasp.a<=, vasp.a>= for absolute comparisons

-- CHANGE: clear stored vasp of an object if a "set" message with an invalid vasp is received

-- CHANGE: cleaned up the process loops (with thread yield) and fixed a corresponding bug in vasp.peaks/valleys

-- CHANGE: warn - but don't stop - on large FFT prime factors

-- FIX: vasp.radd and vasp.rpow had wrong argument initialization

-- FIX: mixfft is not thread-safe... introduced thread mutex for it

-- FIX: vasp.shift: shift argument was modulo frame count... why?

-

-0.0.8:

-- CHANGE: vasp.radio got right outlet for non-radio messages

-- ADD: thread yield for all simple vasp objects

-- CHANGE: vasp.? and vasp.??: made explicit list outlets

-- NEW: vasp.window, vasp.*window, vasp.!window, vasp.*!window, vasp.xwindow, vasp.*xwindow

-- REMOVE: vasp.bevel, vasp.!bevel... obsolete due to new vasp.window objects

-- REMOVE: vasp.noradio (vasp.!radio)... can easily be an abstraction of vasp.radio

-- ADD: minimum length parameter for immediate vasps

-- NEW: vasp.channel? (vasp.c?), vasp.channel (vasp.c) for set/query of buffer channel

-- CHANGE: env - reversed time/value pairs to value/time pairs

-- FIX: vasp.! with length arg outputs its reference on bang

-- FIX: on set message check for invalid vasps

-- FIX: reintroduce test for invalid buffer names (why was it commented out?)

-- FIX: fixe several bugs in vasp.offset= and vasp.frames=

-- CHANGE: allow size 0 (is default) for immediate vasps

-

-0.0.7:

-- CHANGE: vasp.m is now called vasp.multi by default

-- CHANGE: vasp.sync and vasp.multi now have a default argument = 2

-- FIX: vasp.phasor - period length now 1 instead of 2 PI

-- FIX: bug in multi mode (more src, one arg vector)

-- ADD: path specification for help symbols

-- CHANGE: rewrote vasp.imm implementation (still very inefficient!!)

-- ADD: prior message defines differential priority for detached operations

-- ADD: full DFT implementation (radix-2/radix-n) - also support of "to" message

-- ADD: full implementation of vasp.(x)tilt with several interpolation methods

-- CHANGE: vasp.split and vasp.join now have default argument = 2

-- ADD: vasp.offset= (vasp.o=) and vasp.frames= (vasp.f=): search for sample values (eg. zero crossings)

-- ADD: checks and log messages for out-of-range buffer offset and length

-

-0.0.6:

-- reorganized file structure

-- ADD: main names of vasp.n and vasp.n? are now vasp.vector and vasp.vectors?, respectively

-- ADD: main names of vasp.? and vasp.?? are now vasp.list and vasp.nonzero, respectively

-- NEW: vasp.size, vasp.size+, vasp.size? (vasp.s,vasp.s+,vasp.s?) for buffer size manipulation

-- NEW: vasp.peaks? for extremum extraction

-- ADD: shortcuts for vasp.frames,vasp.frames+,vasp.frames? -> vasp.f,vasp.f+,vasp.f?

-- CHANGE: changed vasp.offs,vasp.offs+,vasp.offs? to vasp.offset,vasp.offset+,vasp.offset?

-- ADD: shortcuts for vasp.offset,vasp.offset+,vasp.offset? -> vasp.o,vasp.o+,vasp.o?

-- ADD: envelopes (env ...) as arguments to all operations where vasps are used

-- FIX: default arguments for all binary/anytype operations

-- ADD: vasp.frames* (vasp.f*), vasp.frames/ (vasp.f/), vasp.size* (vasp.s*), vasp.size/ (vasp.s/)

-- ADD: detached operation: operations run as threads, according to detach flag/message

-- CHANGE: vasp.sync has as many outputs as inputs and outputs all input vasps

-- ADD: vasp.! : like vasp but stores the content temporarily (not just the reference)

-- ADD: vasp.copy (vasp.->) and vasp.ccopy (vasp.c->) for instant vasp copying

-- ADD: vasp.radio and vasp.noradio (vasp.!radio) ... filters for radio messages

-- ADD: vasp.fix - bashes NANs to zero, normalizes denormal numbers

-- ADD: double type consisting of 2 additive floats (e.g. "double 1. 1.e-13") for all numeric arguments

-- ADD: vasp.(x)shift - "fill" method/flag defines how to fill shifted areas (0..zero (default),1..none,2..edge value)

-

-0.0.5:

-- FIX: lacking sqrt in [vasp.rmin?],[vasp.rmax?]

-- FIX: (offs >= frames) bug in [vasp.offs?] 

-- ADD: lacking setup of [vasp.!-] and [vasp.c!-] objects 

-- FIX: buggy [vasp.int] code

-- FIX: recognition of integer arguments

-

-0.0.4:

-- CHANGE: vasp.min/max functions so that a vasp length 0 results in 0

-- REMOVED: [vasp.inv], [vasp.cinv].... already replaced by [vasp.!/ 1]

-- FIX: outlet bug in [vasp.?]

-- ADD: right inlet to [vasp]... just like in [float] etc.

-

-0.0.3:

-- restructured the code for future use of break-point lists (aka envelopes) as arguments

-- changed some object's names

-- new objects: vasp.min?, vasp.max? and the likes

-- fixed Max problem with connecting vasp.min,vasp.max right outlet to number boxes

-- fixed right inlet problem for generator and filter objects

-

-0.0.2:

-- vasp.cmin,vasp.cmax - renamed to vasp.rmin,vasp.rmax

-- reversed vasp channel and offset

-- fixed bug with arguments to complex binary operations

-- vasp.chk: included channel check

-- fixed pointer bug in vasp.?

-- fixed pointer increment bug in vasp.fhp

-- changed VecBlock implementation

-

-0.0.1:

-- defined the vasp

-- quick and dirty setup of most functions (non-interruptible)

+Version history:
+
+0.1.3:
+- CHANGE: exchanged #defines for templates in vector calculation loops
+- FIX: bug for binary operations of complex vasps
+
+0.1.2:
+- FIX: bug in vasp.frames* ... wrong argument
+- FIX: another bug in resize of immediate buffers
+- CHANGE: immediates will be initialized to zero, also upon resize to a larger size (by default, later make this optional)
+- ADD: attributes for all basic parameters: ref,to,arg, update,detach,prior (loglvl,argchk,unit later)
+- ADD: vasp.channel,vasp.vector - attribute "index"
+- ADD: vasp.frames,vasp.frames+,vasp.offset,vasp.offset+,vasp.size,vasp.size+,vasp.imm - attribute "frames"
+- ADD: vasp.frames*,vasp.frames/,vasp.size*,vasp.size/ - attribute "factor"
+- ADD: vasp.part - attribute "parts"
+- ADD: vasp.peaks? - attribute "peaks"
+- ADD: vasp.shift,vasp.xshift - attribute "fill"
+- ADD: vasp.tilt,vasp.xtilt - attributes "fill","inter"
+- ADD: vasp.offset=,vasp.frames= - attributes "dir","slope"
+- CHANGE: update m_method_ virtual functions of same objects
+- FIX: immediate vasp garbage collection is now done in separate thread
+- ADD: vasp.size... objects: attribute "keep" determines whether buffer content should be preserved (default = 1)
+- CHANGE: buffer reference checks upon all operations enforced (by new CVasp (checked vasp) class)
+- FIX: immediate vasps can not have frames > size
+- FIX: immediate garbage collection was not triggered in threaded build
+- FIX: changed forgotten standard ASSERT to FLEXT_ASSERT (thanks to Frank Barknecht)
+
+0.1.1:
+- FIX: bug in resize of immediate buffers
+- CHANGE: check (at actual operation!) whether vector exceeds buffer size (GetRVecs, GetCVecs -> corrlen)
+
+0.1.0:
+- NEW: vasp.a<, vasp.a>, vasp.a<=, vasp.a>= for absolute comparisons
+- CHANGE: clear stored vasp of an object if a "set" message with an invalid vasp is received
+- CHANGE: cleaned up the process loops (with thread yield) and fixed a corresponding bug in vasp.peaks/valleys
+- CHANGE: warn - but don't stop - on large FFT prime factors
+- FIX: vasp.radd and vasp.rpow had wrong argument initialization
+- FIX: mixfft is not thread-safe... introduced thread mutex for it
+- FIX: vasp.shift: shift argument was modulo frame count... why?
+
+0.0.8:
+- CHANGE: vasp.radio got right outlet for non-radio messages
+- ADD: thread yield for all simple vasp objects
+- CHANGE: vasp.? and vasp.??: made explicit list outlets
+- NEW: vasp.window, vasp.*window, vasp.!window, vasp.*!window, vasp.xwindow, vasp.*xwindow
+- REMOVE: vasp.bevel, vasp.!bevel... obsolete due to new vasp.window objects
+- REMOVE: vasp.noradio (vasp.!radio)... can easily be an abstraction of vasp.radio
+- ADD: minimum length parameter for immediate vasps
+- NEW: vasp.channel? (vasp.c?), vasp.channel (vasp.c) for set/query of buffer channel
+- CHANGE: env - reversed time/value pairs to value/time pairs
+- FIX: vasp.! with length arg outputs its reference on bang
+- FIX: on set message check for invalid vasps
+- FIX: reintroduce test for invalid buffer names (why was it commented out?)
+- FIX: fixe several bugs in vasp.offset= and vasp.frames=
+- CHANGE: allow size 0 (is default) for immediate vasps
+
+0.0.7:
+- CHANGE: vasp.m is now called vasp.multi by default
+- CHANGE: vasp.sync and vasp.multi now have a default argument = 2
+- FIX: vasp.phasor - period length now 1 instead of 2 PI
+- FIX: bug in multi mode (more src, one arg vector)
+- ADD: path specification for help symbols
+- CHANGE: rewrote vasp.imm implementation (still very inefficient!!)
+- ADD: prior message defines differential priority for detached operations
+- ADD: full DFT implementation (radix-2/radix-n) - also support of "to" message
+- ADD: full implementation of vasp.(x)tilt with several interpolation methods
+- CHANGE: vasp.split and vasp.join now have default argument = 2
+- ADD: vasp.offset= (vasp.o=) and vasp.frames= (vasp.f=): search for sample values (eg. zero crossings)
+- ADD: checks and log messages for out-of-range buffer offset and length
+
+0.0.6:
+- reorganized file structure
+- ADD: main names of vasp.n and vasp.n? are now vasp.vector and vasp.vectors?, respectively
+- ADD: main names of vasp.? and vasp.?? are now vasp.list and vasp.nonzero, respectively
+- NEW: vasp.size, vasp.size+, vasp.size? (vasp.s,vasp.s+,vasp.s?) for buffer size manipulation
+- NEW: vasp.peaks? for extremum extraction
+- ADD: shortcuts for vasp.frames,vasp.frames+,vasp.frames? -> vasp.f,vasp.f+,vasp.f?
+- CHANGE: changed vasp.offs,vasp.offs+,vasp.offs? to vasp.offset,vasp.offset+,vasp.offset?
+- ADD: shortcuts for vasp.offset,vasp.offset+,vasp.offset? -> vasp.o,vasp.o+,vasp.o?
+- ADD: envelopes (env ...) as arguments to all operations where vasps are used
+- FIX: default arguments for all binary/anytype operations
+- ADD: vasp.frames* (vasp.f*), vasp.frames/ (vasp.f/), vasp.size* (vasp.s*), vasp.size/ (vasp.s/)
+- ADD: detached operation: operations run as threads, according to detach flag/message
+- CHANGE: vasp.sync has as many outputs as inputs and outputs all input vasps
+- ADD: vasp.! : like vasp but stores the content temporarily (not just the reference)
+- ADD: vasp.copy (vasp.->) and vasp.ccopy (vasp.c->) for instant vasp copying
+- ADD: vasp.radio and vasp.noradio (vasp.!radio) ... filters for radio messages
+- ADD: vasp.fix - bashes NANs to zero, normalizes denormal numbers
+- ADD: double type consisting of 2 additive floats (e.g. "double 1. 1.e-13") for all numeric arguments
+- ADD: vasp.(x)shift - "fill" method/flag defines how to fill shifted areas (0..zero (default),1..none,2..edge value)
+
+0.0.5:
+- FIX: lacking sqrt in [vasp.rmin?],[vasp.rmax?]
+- FIX: (offs >= frames) bug in [vasp.offs?] 
+- ADD: lacking setup of [vasp.!-] and [vasp.c!-] objects 
+- FIX: buggy [vasp.int] code
+- FIX: recognition of integer arguments
+
+0.0.4:
+- CHANGE: vasp.min/max functions so that a vasp length 0 results in 0
+- REMOVED: [vasp.inv], [vasp.cinv].... already replaced by [vasp.!/ 1]
+- FIX: outlet bug in [vasp.?]
+- ADD: right inlet to [vasp]... just like in [float] etc.
+
+0.0.3:
+- restructured the code for future use of break-point lists (aka envelopes) as arguments
+- changed some object's names
+- new objects: vasp.min?, vasp.max? and the likes
+- fixed Max problem with connecting vasp.min,vasp.max right outlet to number boxes
+- fixed right inlet problem for generator and filter objects
+
+0.0.2:
+- vasp.cmin,vasp.cmax - renamed to vasp.rmin,vasp.rmax
+- reversed vasp channel and offset
+- fixed bug with arguments to complex binary operations
+- vasp.chk: included channel check
+- fixed pointer bug in vasp.?
+- fixed pointer increment bug in vasp.fhp
+- changed VecBlock implementation
+
+0.0.1:
+- defined the vasp
+- quick and dirty setup of most functions (non-interruptible)
diff --git a/externals/grill/vasp/make-files.txt b/externals/grill/vasp/make-files.txt
index e842e8f7..50ecb538 100644
--- a/externals/grill/vasp/make-files.txt
+++ b/externals/grill/vasp/make-files.txt
@@ -8,6 +8,6 @@ main.cpp     obj_offs.cpp    obj_size.cpp   opbase.cpp     ops_cmp.cpp     ops_g
               
 HDRS= \
 arg.h      env.h     opdefs.h   oppermute.h  ops_assign.h  ops_cplx.h     ops_flt.h    ops_resmp.h   ops_wnd.h  vbuffer.h \
-buflib.h   main.h    oploop.h   ops.h        ops_carith.h  ops_dft.h      ops_gen.h    ops_search.h  util.h     vecblk.h \
+buflib.h   main.h    opfuns.h 	oploop.h   ops.h        ops_carith.h  ops_dft.h      ops_gen.h    ops_search.h  util.h     vecblk.h \
 classes.h  opbase.h  opparam.h  ops_arith.h  ops_cmp.h     ops_feature.h  ops_rearr.h  ops_trnsc.h   vasp.h
 
diff --git a/externals/grill/vasp/source/main.cpp b/externals/grill/vasp/source/main.cpp
index 63e71095..504f90ba 100644
--- a/externals/grill/vasp/source/main.cpp
+++ b/externals/grill/vasp/source/main.cpp
@@ -26,7 +26,7 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 
 
 
-const C *VASP_VERSION = "0.1.2";
+const C *VASP_VERSION = "0.1.3pre1";
 
 
 
diff --git a/externals/grill/vasp/source/opbase.cpp b/externals/grill/vasp/source/opbase.cpp
index 989a4952..f3d10ea4 100644
--- a/externals/grill/vasp/source/opbase.cpp
+++ b/externals/grill/vasp/source/opbase.cpp
@@ -84,12 +84,7 @@ Vasp *VaspOp::m_cbin(OpParam &p,CVasp &src,const Argument &arg,CVasp *dst,VecOp:
 	return ret;
 }
 
-
-BL VecOp::_d__run(V fun(S &v,S a),OpParam &p)	{ _D__run(fun,p); }
-BL VecOp::_d__cun(V fun(S &rv,S &iv,S ra,S ia),OpParam &p) { _D__cun(fun,p); }
-BL VecOp::_d__rbin(V fun(S &v,S a,S b),OpParam &p) { _D__rbin(fun,p); }
-BL VecOp::_d__cbin(V fun(S &rv,S &iv,S ra,S ia,S rb,S ib),OpParam &p) { _D__cbin(fun,p); }
-BL VecOp::_d__rop(V fun(S &v,S a,OpParam &p),OpParam &p) { _D__rop(fun,p); }
-BL VecOp::_d__cop(V fun(S &rv,S &iv,S ra,S ia,OpParam &p),OpParam &p) { _D__cop(fun,p); }
-
+#ifdef FLEXT_THREADS
+flext::ThrMutex VecOp::C_base::mtx;
+#endif
 
diff --git a/externals/grill/vasp/source/opbase.h b/externals/grill/vasp/source/opbase.h
index a429b7a6..47cc61cf 100644
--- a/externals/grill/vasp/source/opbase.h
+++ b/externals/grill/vasp/source/opbase.h
@@ -19,12 +19,73 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 namespace VecOp {
 	typedef BL opfun(OpParam &p);
 
-	BL _d__run(V fun(S &v,S a),OpParam &p);
-	BL _d__cun(V fun(S &rv,S &iv,S ra,S ia),OpParam &p);
-	BL _d__rbin(V fun(S &v,S a,S b),OpParam &p);
-	BL _d__cbin(V fun(S &rv,S &iv,S ra,S ia,S rb,S ib),OpParam &p);
-	BL _d__rop(V fun(S &v,S a,OpParam &p),OpParam &p);
-	BL _d__cop(V fun(S &rv,S &iv,S ra,S ia,OpParam &p),OpParam &p);
+    class C_base {
+    public:
+    #ifdef FLEXT_THREADS
+        static flext::ThrMutex mtx;
+        static V Lock() { mtx.Lock(); }
+        static V Unlock() { mtx.Unlock(); }
+    #else
+        static V Lock() {}
+        static V Unlock() {}
+    #endif
+    };
+
+    template<class T> class C_run: public C_base {
+    public: 
+        static BL Do(V f(T &v,T a),OpParam &p) { Lock(); fun = f; _D__run<T,C_run<T> >(p); Unlock(); return true; }
+        static V run(T &v,T a) { fun(v,a); } 
+        static V (*fun)(T &v,T a);
+    };
+    template<class T> V (*C_run<T>::fun)(T &v,T a);
+
+    template<class T> class C_cun: public C_base {
+    public: 
+        static BL Do(V f(T &rv,T &iv,T ra,T ia),OpParam &p) { Lock(); fun = f; _D__cun<T,C_cun<T> >(p); Unlock(); return true; }
+        static V cun(T &rv,T &iv,T ra,T ia) { fun(rv,iv,ra,ia); } 
+        static V (*fun)(T &rv,T &iv,T ra,T ia);
+    };
+    template<class T> V (*C_cun<T>::fun)(T &rv,T &iv,T ra,T ia);
+
+    template<class T> class C_rbin: public C_base {
+    public: 
+        static BL Do(V f(T &v,T a,T b),OpParam &p) { Lock(); fun = f; _D__rbin<T,C_rbin<T> >(p); Unlock(); return true; }
+        static V rbin(T &v,T a,T b) { fun(v,a,b); } 
+        static V (*fun)(T &v,T a,T b);
+    };
+    template<class T> V (*C_rbin<T>::fun)(T &v,T a,T b);
+
+    template<class T> class C_cbin: public C_base {
+    public: 
+        static BL Do(V f(T &rv,T &iv,T ra,T ia,T rb,T ib),OpParam &p) { Lock(); fun = f; _D__cbin<T,C_cbin<T> >(p); Unlock(); return true; }
+        static V cbin(T &rv,T &iv,T ra,T ia,T rb,T ib) { fun(rv,iv,ra,ia,rb,ib); } 
+        static V (*fun)(T &rv,T &iv,T ra,T ia,T rb,T ib);
+    };
+    template<class T> V (*C_cbin<T>::fun)(T &rv,T &iv,T ra,T ia,T rb,T ib);
+
+    template<class T> class C_rop: public C_base {
+    public: 
+        static BL Do(V f(T &v,T a,OpParam &p),OpParam &p) { Lock(); fun = f; _D__rop<T,C_rop<T> >(p); Unlock(); return true; }
+        static V rop(T &v,T a,OpParam &p) { fun(v,a,p); } 
+        static V (*fun)(T &v,T a,OpParam &p);
+    };
+    template<class T> V (*C_rop<T>::fun)(T &v,T a,OpParam &p);
+
+    template<class T> class C_cop: public C_base {
+    public: 
+        static BL Do(V f(T &rv,T &iv,T ra,T ia,OpParam &p),OpParam &p) { Lock(); fun = f; _D__cop<T,C_cop<T> >(p); Unlock(); return true; }
+        static V cop(T &rv,T &iv,T ra,T ia,OpParam &p) { fun(rv,iv,ra,ia,p); } 
+        static V (*fun)(T &rv,T &iv,T ra,T ia,OpParam &p);
+    };
+    template<class T> V (*C_cop<T>::fun)(T &rv,T &iv,T ra,T ia,OpParam &p);
+
+
+    template<class T> BL _d__run(V fun(T &v,T a),OpParam &p)	{ return C_run<T>::Do(fun,p); }
+    template<class T> BL _d__cun(V fun(T &rv,T &iv,T ra,T ia),OpParam &p) { return C_cun<T>::Do(fun,p); }
+    template<class T> BL _d__rbin(V fun(T &v,T a,T b),OpParam &p) { return C_rbin<T>::Do(fun,p); }
+    template<class T> BL _d__cbin(V fun(T &rv,T &iv,T ra,T ia,T rb,T ib),OpParam &p) { return C_cbin<T>::Do(fun,p); }
+    template<class T> BL _d__rop(V fun(T &v,T a,OpParam &p),OpParam &p) { return C_rop<T>::Do(fun,p); }
+    template<class T> BL _d__cop(V fun(T &rv,T &iv,T ra,T ia,OpParam &p),OpParam &p) { return C_cop<T>::Do(fun,p); }
 }
 
 
diff --git a/externals/grill/vasp/source/opdefs.h b/externals/grill/vasp/source/opdefs.h
index bf9101f7..18863f61 100644
--- a/externals/grill/vasp/source/opdefs.h
+++ b/externals/grill/vasp/source/opdefs.h
@@ -12,6 +12,7 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 #define __VASP_OPDEFS_H
 
 #include "oploop.h"
+#include "opbase.h"
 
 #ifdef VASP_CHN1  
 #define _D_ALWAYS1 1
@@ -19,346 +20,375 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 #define _D_ALWAYS1 0
 #endif
 
+namespace VecOp {
+
 /*! \brief skeleton for unary real operations
 */
-#define _D__run(fun,p)											\
-{																\
-	register const S *sr = p.rsdt;								\
-	register S *dr = p.rddt;									\
-	register I i;												\
-	if(sr == dr)												\
-		if(_D_ALWAYS1 || p.rds == 1)							\
-			_D_LOOP(i,p.frames)									\
-			{ fun(*dr,*dr); dr++; }								\
-			_E_LOOP												\
-		else													\
-			_D_LOOP(i,p.frames)									\
-			{ fun(*dr,*dr); dr += p.rds; }						\
-			_E_LOOP												\
-	else														\
-		if(_D_ALWAYS1 || (p.rss == 1 && p.rds == 1))			\
-			_D_LOOP(i,p.frames)									\
-			{ fun(*dr,*sr); sr++,dr++; }						\
-			_E_LOOP												\
-		else													\
-			_D_LOOP(i,p.frames)									\
-			{ fun(*dr,*sr); sr += p.rss,dr += p.rds; }			\
-			_E_LOOP												\
-	return true;												\
+template<class T,V FUN(T &v,T a)> BL _F__run(OpParam &p)										
+{																
+	register const S *sr = p.rsdt;								
+	register S *dr = p.rddt;									
+	register I i;												
+	if(sr == dr)												
+		if(_D_ALWAYS1 || p.rds == 1)							
+			_D_LOOP(i,p.frames)									
+            { FUN(*dr,*dr); dr++; }								
+			_E_LOOP												
+		else													
+			_D_LOOP(i,p.frames)									
+			{ FUN(*dr,*dr); dr += p.rds; }						
+			_E_LOOP												
+	else														
+		if(_D_ALWAYS1 || (p.rss == 1 && p.rds == 1))			
+			_D_LOOP(i,p.frames)									
+			{ FUN(*dr,*sr); sr++,dr++; }						
+			_E_LOOP												
+		else													
+			_D_LOOP(i,p.frames)									
+			{ FUN(*dr,*sr); sr += p.rss,dr += p.rds; }			
+			_E_LOOP												
+	return true;												
 }
 
-#define d__run(fun,p) { return _d__run(fun,p); }
+template<class T,class CL> inline BL _D__run(OpParam &p) { return _F__run<T,CL::run>(p); }
+template<class T,class CL> inline BL d__run(OpParam &p) { return _d__run<T>(CL::run,p); }
 
 /*! \brief skeleton for unary complex operations
 */
-#define _D__cun(fun,p)											\
-{																\
-	register const S *sr = p.rsdt,*si = p.isdt;					\
-	register S *dr = p.rddt,*di = p.iddt;						\
-	register I i;												\
-	if(sr == dr && si == di)									\
-		if(_D_ALWAYS1 || (p.rds == 1 && p.ids == 1))			\
-			_D_LOOP(i,p.frames)									\
-			{ fun(*dr,*di,*dr,*di); dr++,di++; }				\
-			_E_LOOP												\
-		else													\
-			_D_LOOP(i,p.frames)									\
-			{ fun(*dr,*di,*dr,*di); dr += p.rds,di += p.ids; }	\
-			_E_LOOP												\
-	else														\
-		if(_D_ALWAYS1 || (p.rss == 1 && p.iss == 1 && p.rds == 1 && p.ids == 1)) \
-			_D_LOOP(i,p.frames)									\
-			{ fun(*dr,*di,*sr,*si); sr++,si++,dr++,di++; }		\
-			_E_LOOP												\
-		else													\
-			_D_LOOP(i,p.frames)									\
-			{ fun(*dr,*di,*sr,*si);	sr += p.rss,si += p.iss,dr += p.rds,di += p.ids; } \
-			_E_LOOP												\
-	return true;												\
+template<class T,V FUN(T &rv,T &iv,T ra,T ia)> BL _F__cun(OpParam &p)
+{																
+	register const S *sr = p.rsdt,*si = p.isdt;					
+	register S *dr = p.rddt,*di = p.iddt;						
+	register I i;												
+	if(sr == dr && si == di)									
+		if(_D_ALWAYS1 || (p.rds == 1 && p.ids == 1))			
+			_D_LOOP(i,p.frames)									
+			{ FUN(*dr,*di,*dr,*di); dr++,di++; }				
+			_E_LOOP												
+		else													
+			_D_LOOP(i,p.frames)									
+			{ FUN(*dr,*di,*dr,*di); dr += p.rds,di += p.ids; }	
+			_E_LOOP												
+	else														
+		if(_D_ALWAYS1 || (p.rss == 1 && p.iss == 1 && p.rds == 1 && p.ids == 1)) 
+			_D_LOOP(i,p.frames)									
+			{ FUN(*dr,*di,*sr,*si); sr++,si++,dr++,di++; }	
+			_E_LOOP												
+		else													
+			_D_LOOP(i,p.frames)									
+			{ FUN(*dr,*di,*sr,*si);	sr += p.rss,si += p.iss,dr += p.rds,di += p.ids; } 
+			_E_LOOP												
+	return true;												
 }
 
-#define d__cun(fun,p) { return _d__cun(fun,p); }
+template<class T,class CL> inline BL _D__cun(OpParam &p) { return _F__cun<T,CL::cun>(p); }
+template<class T,class CL> inline BL d__cun(OpParam &p) { return _d__cun<T>(CL::cun,p); }
 
 
 /*! \brief skeleton for binary real operations
 */
-#define _D__rbin(fun,p)											\
-{																\
-	register const S *sr = p.rsdt;								\
-	register S *dr = p.rddt;									\
-	register I i;												\
-	if(p.HasArg() && p.arg[0].Is()) {											\
-		switch(p.arg[0].argtp) {									\
-		case OpParam::Arg::arg_v: {								\
-			register const S *ar = p.arg[0].v.rdt;					\
-			if(p.rsdt == p.rddt)									\
-				if(_D_ALWAYS1 || (p.rds == 1 && p.arg[0].v.rs == 1))				\
-					_D_LOOP(i,p.frames)									\
-					{ fun(*dr,*dr,*ar); dr++,ar++; }				\
-					_E_LOOP												\
-				else												\
-					_D_LOOP(i,p.frames)									\
-					{ fun(*dr,*dr,*ar);	dr += p.rds,ar += p.arg[0].v.rs; } \
-					_E_LOOP												\
-			else													\
-				if(_D_ALWAYS1 || (p.rss == 1 && p.rds == 1 && p.arg[0].v.rs == 1))	\
-					_D_LOOP(i,p.frames)									\
-					{ fun(*dr,*sr,*ar); sr++,dr++,ar++; }			\
-					_E_LOOP												\
-				else												\
-					_D_LOOP(i,p.frames)									\
-					{ fun(*dr,*sr,*ar);	sr += p.rss,dr += p.rds,ar += p.arg[0].v.rs; } \
-					_E_LOOP												\
-			break;												\
-		}														\
-		case OpParam::Arg::arg_env: {							\
-			Env::Iter it(*p.arg[0].e.env); it.Init(0);			\
-			if(p.rsdt == p.rddt)									\
-				if(_D_ALWAYS1 || p.rds == 1)										\
-					_D_LOOP(i,p.frames)									\
-					{ fun(*dr,*dr,it.ValFwd(i)); dr++; }			\
-					_E_LOOP												\
-				else												\
-					_D_LOOP(i,p.frames)									\
-					{ fun(*dr,*dr,it.ValFwd(i)); dr += p.rds; }		\
-					_E_LOOP												\
-			else													\
-				if(_D_ALWAYS1 || (p.rss == 1 && p.rds == 1))						\
-					_D_LOOP(i,p.frames)									\
-					{ fun(*dr,*sr,it.ValFwd(i)); sr++,dr++; }		\
-					_E_LOOP												\
-				else												\
-					_D_LOOP(i,p.frames)								 	\
-					{ fun(*dr,*sr,it.ValFwd(i)); sr += p.rss,dr += p.rds; }	\
-					_E_LOOP												\
-			break;												\
-		}														\
-		case OpParam::Arg::arg_x: {							\
-			const R v =  p.arg[0].x.r;							\
-			if(p.rsdt == p.rddt)									\
-				if(_D_ALWAYS1 || p.rds == 1)										\
-					_D_LOOP(i,p.frames)									\
-					{ fun(*dr,*dr,v); dr++; }						\
-					_E_LOOP												\
-				else												\
-					_D_LOOP(i,p.frames)								 	\
-					{ fun(*dr,*dr,v); dr += p.rds; }				\
-					_E_LOOP												\
-			else													\
-				if(_D_ALWAYS1 || (p.rss == 1 && p.rds == 1))						\
-					_D_LOOP(i,p.frames)									\
-					{ fun(*dr,*sr,v); sr++,dr++; }					\
-					_E_LOOP												\
-				else												\
-					_D_LOOP(i,p.frames)								 	\
-					{ fun(*dr,*sr,v); sr += p.rss,dr += p.rds; }	\
-					_E_LOOP												\
-			break;												\
-		}														\
-		}														\
-	}															\
-	else {														\
-		register const S v = p.rbin.arg;						\
-		if(p.rsdt == p.rddt)									\
-			if(_D_ALWAYS1 || p.rds == 1)										\
-				_D_LOOP(i,p.frames)									\
-				{ fun(*dr,*dr,v); dr++; }						\
-				_E_LOOP												\
-			else												\
-				_D_LOOP(i,p.frames)									\
-				{ fun(*dr,*dr,v); dr += p.rds; }				\
-				_E_LOOP												\
-		else													\
-			if(_D_ALWAYS1 || (p.rss == 1 && p.rds == 1))						\
-				_D_LOOP(i,p.frames)									\
-				{ fun(*dr,*sr,v); sr++,dr++; }					\
-				_E_LOOP												\
-			else												\
-				_D_LOOP(i,p.frames)									\
-				{ fun(*dr,*sr,v); sr += p.rss,dr += p.rds; }	\
-				_E_LOOP												\
-	}															\
-	return true;												\
+template<class T,V FUN(T &v,T a,T b)> BL _F__rbin(OpParam &p)					
+{																
+	register const S *sr = p.rsdt;								
+	register S *dr = p.rddt;									
+	register I i;												
+	if(p.HasArg() && p.arg[0].Is()) {							
+		switch(p.arg[0].argtp) {								
+		case OpParam::Arg::arg_v: {								
+			register const S *ar = p.arg[0].v.rdt;				
+			if(p.rsdt == p.rddt)								
+				if(_D_ALWAYS1 || (p.rds == 1 && p.arg[0].v.rs == 1))				
+					_D_LOOP(i,p.frames)									
+                    { FUN(*dr,*dr,*ar); dr++,ar++; }				
+					_E_LOOP												
+				else												
+					_D_LOOP(i,p.frames)									
+					{ FUN(*dr,*dr,*ar);	dr += p.rds,ar += p.arg[0].v.rs; } 
+					_E_LOOP												
+			else													
+				if(_D_ALWAYS1 || (p.rss == 1 && p.rds == 1 && p.arg[0].v.rs == 1))	
+					_D_LOOP(i,p.frames)									
+					{ FUN(*dr,*sr,*ar); sr++,dr++,ar++; }			
+					_E_LOOP												
+				else												
+					_D_LOOP(i,p.frames)									
+					{ FUN(*dr,*sr,*ar);	sr += p.rss,dr += p.rds,ar += p.arg[0].v.rs; } 
+					_E_LOOP												
+			break;												
+		}														
+		case OpParam::Arg::arg_env: {							
+			Env::Iter it(*p.arg[0].e.env); it.Init(0);			
+			if(p.rsdt == p.rddt)								
+				if(_D_ALWAYS1 || p.rds == 1)										
+					_D_LOOP(i,p.frames)									
+					{ FUN(*dr,*dr,it.ValFwd(i)); dr++; }			
+					_E_LOOP												
+				else												
+					_D_LOOP(i,p.frames)									
+					{ FUN(*dr,*dr,it.ValFwd(i)); dr += p.rds; }		
+					_E_LOOP												
+			else													
+				if(_D_ALWAYS1 || (p.rss == 1 && p.rds == 1))						
+					_D_LOOP(i,p.frames)									
+					{ FUN(*dr,*sr,it.ValFwd(i)); sr++,dr++; }		
+					_E_LOOP												
+				else												
+					_D_LOOP(i,p.frames)								 	
+					{ FUN(*dr,*sr,it.ValFwd(i)); sr += p.rss,dr += p.rds; }	
+					_E_LOOP												
+			break;												
+		}														
+		case OpParam::Arg::arg_x: {							
+			const R v =  p.arg[0].x.r;						
+			if(p.rsdt == p.rddt)							
+				if(_D_ALWAYS1 || p.rds == 1)				
+					_D_LOOP(i,p.frames)						
+					{ FUN(*dr,*dr,v); dr++; }			
+					_E_LOOP									
+				else										
+					_D_LOOP(i,p.frames)						
+					{ FUN(*dr,*dr,v); dr += p.rds; }	
+					_E_LOOP									
+			else											
+				if(_D_ALWAYS1 || (p.rss == 1 && p.rds == 1))
+					_D_LOOP(i,p.frames)						
+					{ FUN(*dr,*sr,v); sr++,dr++; }		
+					_E_LOOP									
+				else										
+					_D_LOOP(i,p.frames)						
+                    { FUN(*dr,*sr,v); sr += p.rss,dr += p.rds; }
+					_E_LOOP											
+			break;												
+		}														
+		}														
+	}															
+	else {														
+		register const S v = p.rbin.arg;						
+		if(p.rsdt == p.rddt)									
+			if(_D_ALWAYS1 || p.rds == 1)						
+				_D_LOOP(i,p.frames)									
+				{ FUN(*dr,*dr,v); dr++; }						
+				_E_LOOP												
+			else												
+				_D_LOOP(i,p.frames)									
+				{ FUN(*dr,*dr,v); dr += p.rds; }				
+				_E_LOOP												
+		else													
+			if(_D_ALWAYS1 || (p.rss == 1 && p.rds == 1))		
+				_D_LOOP(i,p.frames)									
+				{ FUN(*dr,*sr,v); sr++,dr++; }					
+				_E_LOOP												
+			else												
+				_D_LOOP(i,p.frames)									
+				{ FUN(*dr,*sr,v); sr += p.rss,dr += p.rds; }	
+				_E_LOOP												
+	}															
+	return true;												
 }
 
-#define d__rbin(fun,p) { return _d__rbin(fun,p); }
+template<class T,class CL> inline BL _D__rbin(OpParam &p) { return _F__rbin<T,CL::rbin>(p); }
+template<class T,class CL> inline BL d__rbin(OpParam &p) { return _d__rbin<T>(CL::rbin,p); }
 
 
 /*! \brief skeleton for binary complex operations
 */
-#define _D__cbin(fun,p)											\
-{																\
-	register const S *sr = p.rsdt,*si = p.isdt;					\
-	register S *dr = p.rddt,*di = p.iddt;						\
-	register I i;												\
-	if(p.HasArg() && p.arg[0].Is()) {											\
-		switch(p.arg[0].argtp) {									\
-		case OpParam::Arg::arg_v: {									\
-			register const S *ar = p.arg[0].v.rdt,*ai = p.arg[0].v.idt;				\
-			if(ai)													\
-				if(sr == dr && si == di)							\
-					if(_D_ALWAYS1 || (p.rds == 1 && p.ids == 1 && p.arg[0].v.rs == 1 && p.arg[0].v.is == 1)) \
-						_D_LOOP(i,p.frames)				 \
-						{ fun(*dr,*di,*dr,*di,*ar,*ai);	dr++,di++,ar++,ai++; }		\
-						_E_LOOP												\
-					else											\
-						_D_LOOP(i,p.frames)				 \
-						{ fun(*dr,*di,*dr,*di,*ar,*ai);	dr += p.rds,di += p.ids,ar += p.arg[0].v.rs,ai += p.arg[0].v.is; }		\
-						_E_LOOP												\
-				else												\
-					_D_LOOP(i,p.frames)				 \
-					{ fun(*dr,*di,*sr,*si,*ar,*ai);	sr += p.rss,si += p.iss,dr += p.rds,di += p.ids,ar += p.arg[0].v.rs,ai += p.arg[0].v.is; }			\
-					_E_LOOP												\
-			else													\
-				if(sr == dr && si == di)							\
-					_D_LOOP(i,p.frames)				 \
-					{ fun(*dr,*di,*dr,*di,*ar,0); dr += p.rds,di += p.ids,ar += p.arg[0].v.rs; }	\
-					_E_LOOP												\
-				else												\
-					_D_LOOP(i,p.frames)				 \
-					{ fun(*dr,*di,*sr,*si,*ar,0); sr += p.rss,si += p.iss,dr += p.rds,di += p.ids,ar += p.arg[0].v.rs; } \
-					_E_LOOP												\
-			break;												\
-		}														\
-		case OpParam::Arg::arg_env: {									\
-			Env::Iter it(*p.arg[0].e.env); it.Init(0);				\
-			if(sr == dr && si == di)							\
-				if(_D_ALWAYS1 || (p.rds == 1 && p.ids == 1)) \
-					_D_LOOP(i,p.frames)				 \
-					{ fun(*dr,*di,*dr,*di,it.ValFwd(i),0); dr++,di++; }			\
-					_E_LOOP												\
-				else											\
-					_D_LOOP(i,p.frames)				 \
-					{ fun(*dr,*di,*dr,*di,it.ValFwd(i),0); dr += p.rds,di += p.ids; }			\
-					_E_LOOP												\
-			else												\
-				_D_LOOP(i,p.frames)				 \
-					{ fun(*dr,*di,*sr,*si,it.ValFwd(i),0); sr += p.rss,si += p.iss,dr += p.rds,di += p.ids; }	\
-				_E_LOOP												\
-			break;												\
-		}														\
-		case OpParam::Arg::arg_x: {									\
-			register const R ar = p.arg[0].x.r,ai = p.arg[0].x.i;				\
-			if(sr == dr && si == di)							\
-				if(_D_ALWAYS1 || (p.rds == 1 && p.ids == 1)) \
-					_D_LOOP(i,p.frames)				 \
-					{ fun(*dr,*di,*dr,*di,ar,ai); dr++,di++; }			\
-					_E_LOOP												\
-				else											\
-					_D_LOOP(i,p.frames)				 \
-					{ fun(*dr,*di,*dr,*di,ar,ai); dr += p.rds,di += p.ids; }			\
-					_E_LOOP												\
-			else												\
-				_D_LOOP(i,p.frames)				 \
-				{ fun(*dr,*di,*sr,*si,ar,ai); sr += p.rss,si += p.iss,dr += p.rds,di += p.ids; }	\
-				_E_LOOP												\
-			break;												\
-		}														\
-		}														\
-	}															\
-	else {														\
-		register const S rv = p.cbin.rarg,iv = p.cbin.iarg;		\
-		if(sr == dr && si == di)								\
-			if(_D_ALWAYS1 || (p.rds == 1 && p.ids == 1))						\
-				_D_LOOP(i,p.frames)						\
-				{ fun(*dr,*di,*dr,*di,rv,iv); dr++,di++; }					\
-				_E_LOOP												\
-			else												\
-				_D_LOOP(i,p.frames)				 \
-				{ fun(*dr,*di,*dr,*di,rv,iv); dr += p.rds,di += p.ids; }					\
-				_E_LOOP												\
-		else													\
-			if(_D_ALWAYS1 || (p.rds == 1 && p.ids == 1 && p.rss == 1 && p.iss == 1)) \
-				_D_LOOP(i,p.frames)				 \
-				{ fun(*dr,*di,*sr,*si,rv,iv); sr++,si++,dr++,di++; }					\
-				_E_LOOP												\
-			else												\
-				_D_LOOP(i,p.frames)				 \
-				{ fun(*dr,*di,*sr,*si,rv,iv); sr += p.rss,si += p.iss,dr += p.rds,di += p.ids; }					\
-				_E_LOOP												\
-	}															\
-	return true;												\
+template<class T,V FUN(T &rv,T &iv,T ra,T ia,T rb,T ib)> BL _F__cbin(OpParam &p)							
+{																
+	register const S *sr = p.rsdt,*si = p.isdt;					
+	register S *dr = p.rddt,*di = p.iddt;						
+	register I i;												
+	if(p.HasArg() && p.arg[0].Is()) {							
+		switch(p.arg[0].argtp) {									
+		case OpParam::Arg::arg_v: {									
+			register const S *ar = p.arg[0].v.rdt,*ai = p.arg[0].v.idt;				
+			if(ai)													
+				if(sr == dr && si == di)							
+					if(_D_ALWAYS1 || (p.rds == 1 && p.ids == 1 && p.arg[0].v.rs == 1 && p.arg[0].v.is == 1)) 
+						_D_LOOP(i,p.frames)				 
+						{ FUN(*dr,*di,*dr,*di,*ar,*ai);	dr++,di++,ar++,ai++; }		
+						_E_LOOP												
+					else											
+						_D_LOOP(i,p.frames)				 
+						{ FUN(*dr,*di,*dr,*di,*ar,*ai);	dr += p.rds,di += p.ids,ar += p.arg[0].v.rs,ai += p.arg[0].v.is; }		
+						_E_LOOP												
+				else												
+					_D_LOOP(i,p.frames)				 
+					{ FUN(*dr,*di,*sr,*si,*ar,*ai);	sr += p.rss,si += p.iss,dr += p.rds,di += p.ids,ar += p.arg[0].v.rs,ai += p.arg[0].v.is; }			
+					_E_LOOP												
+			else													
+				if(sr == dr && si == di)							
+					_D_LOOP(i,p.frames)				 
+					{ FUN(*dr,*di,*dr,*di,*ar,0); dr += p.rds,di += p.ids,ar += p.arg[0].v.rs; }	
+					_E_LOOP												
+				else												
+					_D_LOOP(i,p.frames)				 
+					{ FUN(*dr,*di,*sr,*si,*ar,0); sr += p.rss,si += p.iss,dr += p.rds,di += p.ids,ar += p.arg[0].v.rs; } 
+					_E_LOOP												
+			break;												
+		}														
+		case OpParam::Arg::arg_env: {							
+			Env::Iter it(*p.arg[0].e.env); it.Init(0);			
+			if(sr == dr && si == di)							
+				if(_D_ALWAYS1 || (p.rds == 1 && p.ids == 1)) 
+					_D_LOOP(i,p.frames)				 
+					{ FUN(*dr,*di,*dr,*di,it.ValFwd(i),0); dr++,di++; }			
+					_E_LOOP												
+				else											
+					_D_LOOP(i,p.frames)				 
+					{ FUN(*dr,*di,*dr,*di,it.ValFwd(i),0); dr += p.rds,di += p.ids; }			
+					_E_LOOP												
+			else												
+				_D_LOOP(i,p.frames)				 
+					{ FUN(*dr,*di,*sr,*si,it.ValFwd(i),0); sr += p.rss,si += p.iss,dr += p.rds,di += p.ids; }	
+				_E_LOOP												
+			break;												
+		}														
+		case OpParam::Arg::arg_x: {								
+			register const R ar = p.arg[0].x.r,ai = p.arg[0].x.i;				
+			if(sr == dr && si == di)							
+				if(_D_ALWAYS1 || (p.rds == 1 && p.ids == 1)) 
+					_D_LOOP(i,p.frames)				 
+					{ FUN(*dr,*di,*dr,*di,ar,ai); dr++,di++; }			
+					_E_LOOP												
+				else											
+					_D_LOOP(i,p.frames)				 
+					{ FUN(*dr,*di,*dr,*di,ar,ai); dr += p.rds,di += p.ids; }			
+					_E_LOOP												
+			else												
+				_D_LOOP(i,p.frames)				 
+				{ FUN(*dr,*di,*sr,*si,ar,ai); sr += p.rss,si += p.iss,dr += p.rds,di += p.ids; }	
+				_E_LOOP												
+			break;												
+		}														
+		}														
+	}															
+	else {														
+		register const S rv = p.cbin.rarg,iv = p.cbin.iarg;		
+		if(sr == dr && si == di)								
+			if(_D_ALWAYS1 || (p.rds == 1 && p.ids == 1))		
+				_D_LOOP(i,p.frames)						
+				{ FUN(*dr,*di,*dr,*di,rv,iv); dr++,di++; }					
+				_E_LOOP												
+			else												
+				_D_LOOP(i,p.frames)				 
+				{ FUN(*dr,*di,*dr,*di,rv,iv); dr += p.rds,di += p.ids; }					
+				_E_LOOP												
+		else													
+			if(_D_ALWAYS1 || (p.rds == 1 && p.ids == 1 && p.rss == 1 && p.iss == 1)) 
+				_D_LOOP(i,p.frames)				 
+				{ FUN(*dr,*di,*sr,*si,rv,iv); sr++,si++,dr++,di++; }					
+				_E_LOOP												
+			else												
+				_D_LOOP(i,p.frames)				 
+				{ FUN(*dr,*di,*sr,*si,rv,iv); sr += p.rss,si += p.iss,dr += p.rds,di += p.ids; }					
+				_E_LOOP												
+	}															
+	return true;												
 }
 
-#define d__cbin(fun,p) { return _d__cbin(fun,p); }
+template<class T,class CL> inline BL _D__cbin(OpParam &p) { return _F__cbin<T,CL::cbin>(p); }			
+template<class T,class CL> inline BL d__cbin(OpParam &p) { return _d__cbin<T>(CL::cbin,p); }
 
 
 /*! \brief skeleton for real operations with parameter block
 */
-#define _D__rop(fun,p)											\
-{																\
-	register const S *sr = p.rsdt;								\
-	register S *dr = p.rddt;									\
-	register I i;												\
-	if(sr == dr)												\
-		if(_D_ALWAYS1 || p.rds == 1)											\
-			_D_LOOP(i,p.frames)										\
-			{ fun(*dr,*dr,p); dr++; }							\
-			_E_LOOP												\
-		else													\
-			_D_LOOP(i,p.frames)										\
-			{ fun(*dr,*dr,p); dr += p.rds; }					\
-			_E_LOOP												\
-	else														\
-		if(_D_ALWAYS1 || (p.rss == 1 && p.rds == 1))							\
-			_D_LOOP(i,p.frames)										\
-			{ fun(*dr,*sr,p); sr++,dr++; }						\
-			_E_LOOP												\
-		else													\
-			_D_LOOP(i,p.frames)										\
-			{ fun(*dr,*sr,p); sr += p.rss,dr += p.rds; }		\
-			_E_LOOP												\
-	return true;												\
+template<class T,V FUN(T &v,T r,OpParam &p)> BL _F__rop(OpParam &p)						
+{																
+	register const S *sr = p.rsdt;								
+	register S *dr = p.rddt;									
+	register I i;												
+	if(sr == dr)												
+		if(_D_ALWAYS1 || p.rds == 1)											
+			_D_LOOP(i,p.frames)										
+            { FUN(*dr,*dr,p); dr++; }							
+			_E_LOOP												
+		else													
+			_D_LOOP(i,p.frames)									
+			{ FUN(*dr,*dr,p); dr += p.rds; }					
+			_E_LOOP												
+	else														
+		if(_D_ALWAYS1 || (p.rss == 1 && p.rds == 1))			
+			_D_LOOP(i,p.frames)									
+			{ FUN(*dr,*sr,p); sr++,dr++; }					
+			_E_LOOP												
+		else													
+			_D_LOOP(i,p.frames)									
+			{ FUN(*dr,*sr,p); sr += p.rss,dr += p.rds; }		
+			_E_LOOP												
+	return true;												
 }
 
-#define d__rop(fun,p) { return _d__rop(fun,p); }
+template<class T,class CL> inline BL _D__rop(OpParam &p) { return _F__rop<T,CL::rop>(p); }			
+template<class T,class CL> inline BL d__rop(OpParam &p) { return _d__rop<T>(CL::rop,p); }
 
 
 /*! \brief skeleton for complex operations with parameter block
 */
-#define _D__cop(fun,p)											\
-{																\
-	register const S *sr = p.rsdt,*si = p.isdt;					\
-	register S *dr = p.rddt,*di = p.iddt;						\
-	register I i;												\
-	if(sr == dr && si == di)									\
-		if(_D_ALWAYS1 || (p.rds == 1 && p.ids == 1))							\
-			_D_LOOP(i,p.frames)										\
-			{ fun(*dr,*di,*dr,*di,p); dr++,di++; }				\
-			_E_LOOP												\
-		else													\
-			_D_LOOP(i,p.frames)										\
-			{ fun(*dr,*di,*dr,*di,p); dr += p.rds,di += p.ids; } \
-			_E_LOOP												\
-	else														\
-		if(_D_ALWAYS1 || (p.rss == 1 && p.iss == 1 && p.rds == 1 && p.ids == 1)) \
-			_D_LOOP(i,p.frames)				 \
-			{ fun(*dr,*di,*sr,*si,p); sr++,si++,dr++,di++; }	\
-			_E_LOOP												\
-		else													\
-			_D_LOOP(i,p.frames)				 \
-			{ fun(*dr,*di,*sr,*si,p); sr += p.rss,si += p.iss,dr += p.rds,di += p.ids; }	\
-			_E_LOOP												\
-	return true;												\
+template<class T,V FUN(T &rv,T &iv,T ra,T ia,OpParam &p)> BL _F__cop(OpParam &p)						
+{																
+	register const S *sr = p.rsdt,*si = p.isdt;					
+	register S *dr = p.rddt,*di = p.iddt;						
+	register I i;												
+	if(sr == dr && si == di)									
+		if(_D_ALWAYS1 || (p.rds == 1 && p.ids == 1))			
+			_D_LOOP(i,p.frames)									
+			{ FUN(*dr,*di,*dr,*di,p); dr++,di++; }			
+			_E_LOOP												
+		else													
+			_D_LOOP(i,p.frames)									
+			{ FUN(*dr,*di,*dr,*di,p); dr += p.rds,di += p.ids; } 
+			_E_LOOP												
+	else														
+		if(_D_ALWAYS1 || (p.rss == 1 && p.iss == 1 && p.rds == 1 && p.ids == 1)) 
+			_D_LOOP(i,p.frames)				 
+			{ FUN(*dr,*di,*sr,*si,p); sr++,si++,dr++,di++; }	
+			_E_LOOP												
+		else													
+			_D_LOOP(i,p.frames)				 
+			{ FUN(*dr,*di,*sr,*si,p); sr += p.rss,si += p.iss,dr += p.rds,di += p.ids; }	
+			_E_LOOP												
+	return true;												
 }
 
-#define d__cop(fun,p) { return _d__cop(fun,p); }
+template<class T,class CL> inline BL _D__cop(OpParam &p) { return _F__cop<T,CL::cop>(p); }
+template<class T,class CL> inline BL d__cop(OpParam &p) { return _d__cop<T>(CL::cop,p); }
+
+template<class T,V FUN(T &v,T a)> inline BL f__run(OpParam &p) { return _d__run(FUN,p); }
+template<class T,V FUN(T &rv,T &iv,T ra,T ia)> inline BL f__cun(OpParam &p) { return _d__cun(FUN,p); }
+template<class T,V FUN(T &v,T a,T b)> inline BL f__rbin(OpParam &p) { return _d__rbin(FUN,p); }
+template<class T,V FUN(T &rv,T &iv,T ra,T ia,T rb,T ib)> inline BL f__cbin(OpParam &p) { return _d__cbin(FUN,p); }
+template<class T,V FUN(T &v,T r,OpParam &p)> inline BL f__rop(OpParam &p) { return _d__rop(FUN,p); }
+template<class T,V FUN(T &rv,T &iv,T ra,T ia,OpParam &p)> inline BL f__cop(OpParam &p) { return _d__cop(FUN,p); }
 
 
 #ifdef VASP_COMPACT
-#define D__run(fun,p) d__run(fun,p)
-#define D__cun(fun,p) d__cun(fun,p)
-#define D__rbin(fun,p) d__rbin(fun,p)
-#define D__cbin(fun,p) d__cbin(fun,p)
-#define D__rop(fun,p) d__rop(fun,p)
-#define D__cop(fun,p) d__cop(fun,p)
+template<class T,class CL> inline BL D__run(OpParam &p) { return _d__run<T>(CL::run,p); }
+template<class T,class CL> inline BL D__cun(OpParam &p) { return _d__cun<T>(CL::cun,p); }
+template<class T,class CL> inline BL D__rbin(OpParam &p) { return _d__rbin<T>(CL::rbin,p); }
+template<class T,class CL> inline BL D__cbin(OpParam &p) { return _d__cbin<T>(CL::cbin,p); }
+template<class T,class CL> inline BL D__rop(OpParam &p) { return _d__rop<T>(CL::rop,p); }
+template<class T,class CL> inline BL D__cop(OpParam &p) { return _d__cop<T>(CL::cop,p); }
+template<class T,V FUN(T &v,T a)> inline BL F__run(OpParam &p) { return _d__run(FUN,p); }
+template<class T,V FUN(T &rv,T &iv,T ra,T ia)> inline BL F__cun(OpParam &p) { return _d__cun<T>(FUN,p); }
+template<class T,V FUN(T &v,T a,T b)> inline BL F__rbin(OpParam &p) { return _d__rbin<T>(FUN,p); }
+template<class T,V FUN(T &rv,T &iv,T ra,T ia,T rb,T ib)> inline BL F__cbin(OpParam &p) { return _d__cbin<T>(FUN,p); }
+template<class T,V FUN(T &v,T r,OpParam &p)> inline BL F__rop(OpParam &p) { return _d__rop<T>(FUN,p); }
+template<class T,V FUN(T &rv,T &iv,T ra,T ia,OpParam &p)> inline BL F__cop(OpParam &p) { return _d__cop<T>(FUN,p); }
 #else
-#define D__run(fun,p) _D__run(fun,p)
-#define D__cun(fun,p) _D__cun(fun,p)
-#define D__rbin(fun,p) _D__rbin(fun,p)
-#define D__cbin(fun,p) _D__cbin(fun,p)
-#define D__rop(fun,p) _D__rop(fun,p)
-#define D__cop(fun,p) _D__cop(fun,p)
+template<class T,class CL> inline BL D__run(OpParam &p) { return _D__run<T,CL>(p); }
+template<class T,class CL> inline BL D__cun(OpParam &p) { return _D__cun<T,CL>(p); }
+template<class T,class CL> inline BL D__rbin(OpParam &p) { return _D__rbin<T,CL>(p); }
+template<class T,class CL> inline BL D__cbin(OpParam &p) { return _D__cbin<T,CL>(p); }
+template<class T,class CL> inline BL D__rop(OpParam &p) { return _D__rop<T,CL>(p); }
+template<class T,class CL> inline BL D__cop(OpParam &p) { return _D__cop<T,CL>(p); }
+template<class T,V FUN(T &v,T a)> inline BL F__run(OpParam &p) { return _F__run<T,FUN>(p); }
+template<class T,V FUN(T &rv,T &iv,T ra,T ia)> inline BL F__cun(OpParam &p) { return _F__cun<T,FUN>(p); }
+template<class T,V FUN(T &v,T a,T b)> inline BL F__rbin(OpParam &p) { return _F__rbin<T,FUN>(p); }
+template<class T,V FUN(T &rv,T &iv,T ra,T ia,T rb,T ib)> inline BL F__cbin(OpParam &p) { return _F__cbin<T,FUN>(p); }
+template<class T,V FUN(T &v,T r,OpParam &p)> inline BL F__rop(OpParam &p) { return _F__rop<T,FUN>(p); }
+template<class T,V FUN(T &rv,T &iv,T ra,T ia,OpParam &p)> inline BL F__cop(OpParam &p) { return _F__cop<T,FUN>(p); }
 #endif
 
+} // namespace VecOp
+
 #endif
diff --git a/externals/grill/vasp/source/opfuns.h b/externals/grill/vasp/source/opfuns.h
new file mode 100644
index 00000000..99fca899
--- /dev/null
+++ b/externals/grill/vasp/source/opfuns.h
@@ -0,0 +1,378 @@
+/* 
+
+VASP modular - vector assembling signal processor / objects for Max/MSP and PD
+
+Copyright (c) 2002 Thomas Grill (xovo@gmx.net)
+For information on usage and redistribution, and for a DISCLAIMER OF ALL
+WARRANTIES, see the file, "license.txt," in this distribution.  
+
+*/
+
+#ifndef __VASP_OPFUNS_H
+#define __VASP_OPFUNS_H
+
+#include "opdefs.h"
+
+namespace VecOp {
+
+    // assignment
+
+    template<class T> class f_copy { 
+    public: 
+        static V run(T &v,T a) { v = a; }
+        static V cun(T &rv,T &iv,T ra,T ia) { rv = ra,iv = ia; } 
+    };
+
+    template<class T> class f_set { 
+    public: 
+        static V rbin(T &v,T,T b) { v = b; }
+        static V cbin(T &rv,T &iv,T,T,T rb,T ib) { rv = rb,iv = ib; } 
+    };
+
+    // arithmetic
+
+    template<class T> class f_add {
+    public: 
+        static V rbin(T &v,T a,T b) { v = a+b; }
+        static V cbin(T &rv,T &iv,T ra,T ia,T rb,T ib) { rv = ra+rb,iv = ia+ib; }
+    };
+
+    template<class T> class f_sub {
+    public: 
+        static V rbin(T &v,T a,T b) { v = a-b; }
+        static V cbin(T &rv,T &iv,T ra,T ia,T rb,T ib) { rv = ra-rb,iv = ia-ib; }
+    };
+
+    template<class T> class f_subr {
+    public: 
+        static V rbin(T &v,T a,T b) { v = b-a; }
+        static V cbin(T &rv,T &iv,T ra,T ia,T rb,T ib) { rv = rb-ra,iv = ib-ia; }
+    };
+
+    template<class T> class f_mul {
+    public: 
+        static V rbin(T &v,T a,T b) { v = a*b; }
+        static V cbin(T &rv,T &iv,T ra,T ia,T rb,T ib) { rv = ra*rb-ia*ib, iv = ra*ib+rb*ia; }
+    };
+
+    template<class T> class f_div {
+    public: 
+        static V rbin(T &v,T a,T b) { v = a/b; }
+
+        static V cbin(T &rv,T &iv,T ra,T ia,T rb,T ib) 
+        { 
+	        register const R den = sqabs(rb,ib);
+	        rv = (ra*rb+ia*ib)/den;
+	        iv = (ia*rb-ra*ib)/den;
+        }
+    };
+
+    template<class T> class f_divr {
+    public: 
+        static V rbin(T &v,T a,T b) { v = b/a; }
+
+        static V cbin(T &rv,T &iv,T ra,T ia,T rb,T ib)
+        { 
+	        register const R den = sqabs(ra,ia);
+	        rv = (rb*ra+ib*ia)/den;
+	        iv = (ib*ra-rb*ia)/den;
+        }
+    };
+
+    template<class T> class f_mod { 
+    public: 
+        static V rbin(T &v,T a,T b) { v = fmod(a,b); } 
+    };
+
+    template<class T> class f_abs {
+    public: 
+        static V run(T &v,T a) { v = fabs(a); }
+        static V cun(T &rv,T &iv,T ra,T ia) { rv = sqrt(ra*ra+ia*ia),iv = 0; }
+    };
+
+    template<class T> class f_sign { 
+    public: 
+        static V run(T &v,T a) { v = (a == 0?0:(a < 0?-1.:1.)); } 
+    };
+
+    template<class T> class f_sqr {
+    public: 
+        static V run(T &v,T a) { v = a*a; } 
+        static V cun(T &rv,T &iv,T ra,T ia) { rv = ra*ra-ia*ia; iv = ra*ia*2; }
+    };
+
+    template<class T> class f_ssqr { 
+    public: 
+        static V run(T &v,T a) { v = a*fabs(a); } 
+    };
+
+
+    // transcendent
+
+    template<class T> class f_powi {
+    public: 
+        static V cop(T &rv,T &iv,T ra,T ia,OpParam &p) 
+        { 
+	        register const I powi = p.ibin.arg;
+            register S rt,it; f_sqr<T>::cun(rt,it,ra,ia);
+            for(I i = 2; i < powi; ++i) f_mul<T>::cbin(rt,it,rt,it,ra,ia);
+	        rv = rt,iv = it;
+        } 
+    };
+
+    template<class T> class f_pow {
+    public: 
+        static V rbin(T &v,T a,T b) { v = pow(fabs(a),b)*sgn(a); } 
+
+        static V cbin(T &rv,T &iv,T ra,T ia,T rb,T) 
+        { 
+	        register const R _abs = sqrt(sqabs(ra,ia));
+	        if(_abs) {
+		        register const R _p = pow(_abs,rb)/_abs;
+		        rv = _p*ra,iv = _p*ia;
+	        }
+	        else
+		        rv = iv = 0;
+        } 
+    };
+
+    template<class T> class f_sqrt {
+    public: 
+        static V run(T &v,T a) { v = sqrt(fabs(a)); } 
+    };
+
+    template<class T> class f_ssqrt {
+    public:
+        static V run(T &v,T a) { v = sqrt(fabs(a))*sgn(a); } 
+    };
+
+
+    template<class T> class f_exp {
+    public: 
+        static V run(T &v,T a) { v = exp(a); } 
+    };
+
+    template<class T> class f_log {
+    public: 
+        static V run(T &v,T a) { v = log(a); }  // \todo detect NANs
+    };
+
+    // comparisons
+
+    template<class T> class f_lwr {
+    public: 
+        static V rbin(T &v,T a,T b) { v = a < b?1:0; }
+    };
+
+    template<class T> class f_gtr {
+    public: 
+        static V rbin(T &v,T a,T b) { v = a > b?1:0; }
+    };
+
+    template<class T> class f_alwr {
+    public: 
+        static V rbin(T &v,T a,T b) { v = fabs(a) < fabs(b)?1:0; }
+    };
+
+    template<class T> class f_agtr {
+    public: 
+        static V rbin(T &v,T a,T b) { v = fabs(a) > fabs(b)?1:0; }
+    };
+
+    template<class T> class f_leq {
+    public: 
+        static V rbin(T &v,T a,T b) { v = a <= b?1:0; }
+    };
+
+    template<class T> class f_geq {
+    public: 
+        static V rbin(T &v,T a,T b) { v = a >= b?1:0; }
+    };
+
+    template<class T> class f_aleq {
+    public: 
+        static V rbin(T &v,T a,T b) { v = fabs(a) <= fabs(b)?1:0; }
+    };
+
+    template<class T> class f_ageq {
+    public: 
+        static V rbin(T &v,T a,T b) { v = fabs(a) >= fabs(b)?1:0; }
+    };
+
+    template<class T> class f_equ {
+    public: 
+        static V rbin(T &v,T a,T b) { v = a == b?1:0; }
+    };
+
+    template<class T> class f_neq {
+    public: 
+        static V rbin(T &v,T a,T b) { v = a != b?1:0; }
+    };
+
+    // min/max
+
+    template<class T> class f_min {
+    public: 
+        static V rbin(T &v,T a,T b) { v = a < b?a:b; }
+
+        static V cbin(T &rv,T &iv,T ra,T ia,T rb,T ib) 
+        { 
+	        if(sqabs(ra,ia) < sqabs(rb,ib))	rv = ra,iv = ia; 
+	        else rv = rb,iv = ib; 
+        }
+    };
+
+    template<class T> class f_max {
+    public: 
+        static V rbin(T &v,T a,T b) { v = a > b?a:b; }
+
+        static V cbin(T &rv,T &iv,T ra,T ia,T rb,T ib) 
+        { 
+	        if(sqabs(ra,ia) > sqabs(rb,ib))	rv = ra,iv = ia; 
+	        else rv = rb,iv = ib; 
+        }
+    };
+
+    template<class T> class f_minmax {
+    public:
+        static V cun(T &rv,T &iv,T ra,T ia) 
+        { 
+	        if(ra < ia)	rv = ra,iv = ia; 
+	        else rv = ia,iv = ra; 
+        } 
+    };
+
+    template<class T> class f_minq {
+    public: 
+        static V rop(T &,T ra,OpParam &p) 
+        { 
+	        if(ra < p.norm.minmax) p.norm.minmax = ra; 
+        } 
+
+        static V cop(T &,T &,T ra,T ia,OpParam &p) 
+        { 
+	        register T s = sqabs(ra,ia); 
+	        if(s < p.norm.minmax) p.norm.minmax = s; 
+        } 
+    };
+
+    template<class T> class f_maxq {
+    public: 
+        static V rop(T &,T ra,OpParam &p) 
+        { 
+	        if(ra > p.norm.minmax) p.norm.minmax = ra; 
+        } 
+
+        static V cop(T &,T &,T ra,T ia,OpParam &p) 
+        { 
+	        register T s = sqabs(ra,ia); 
+	        if(s > p.norm.minmax) p.norm.minmax = s; 
+        } 
+    };
+
+    template<class T> V f_aminq(T &,T ra,OpParam &p) 
+    { 
+	    register T s = fabs(ra); 
+	    if(s < p.norm.minmax) p.norm.minmax = s; 
+    } 
+
+    template<class T> V f_amaxq(T &,T ra,OpParam &p) 
+    { 
+	    register T s = fabs(ra); 
+	    if(s > p.norm.minmax) p.norm.minmax = s; 
+    } 
+
+
+    // gating
+
+    template<class T> class f_gate {
+    public:
+        static V rbin(T &rv,T ra,T rb) { rv = fabs(ra) >= rb?ra:0; } 
+
+        static V cbin(T &rv,T &iv,T ra,T ia,T rb,T) 
+        { 
+	        register const T _abs = sqabs(ra,ia);
+
+	        if(_abs >= rb*rb) rv = ra,iv = ia;
+	        else rv = iv = 0;
+        } 
+    };
+
+    template<class T> class f_igate {
+    public:
+        static V rbin(T &rv,T ra,T rb) { rv = fabs(ra) <= rb?ra:0; } 
+
+        static V cbin(T &rv,T &iv,T ra,T ia,T rb,T) 
+        { 
+	        register const T _abs = sqabs(ra,ia);
+
+	        if(_abs <= rb*rb) rv = ra,iv = ia;
+	        else rv = iv = 0;
+        } 
+    };
+    
+    // complex
+
+    template<class T> class f_norm {
+    public:
+        static V cun(T &rv,T &iv,T ra,T ia) 
+        { 
+	        register T f = sqabs(ra,ia);
+	        if(f) { f = 1./sqrt(f); rv = ra*f,iv = ia*f; }
+	        else rv = iv = 0;
+        }
+    };
+
+    template<class T> class f_conj {
+    public:
+        static V cun(T &,T &iv,T,T ia) { iv = -ia; }
+    };
+
+    template<class T> class f_polar {
+    public:
+        static V cun(T &rv,T &iv,T ra,T ia) { rv = sqrt(sqabs(ra,ia)),iv = arg(ra,ia); }
+    };
+
+    template<class T> class f_rect {
+    public:
+        static V cun(T &rv,T &iv,T ra,T ia) { rv = ra*cos(ia),iv = ra*sin(ia); }
+    };
+
+    template<class T> class f_radd {
+    public:
+        static V cbin(T &rv,T &iv,T ra,T ia,T rb,T) 
+        { 
+	        register const R _abs = sqrt(sqabs(ra,ia))+rb;
+	        register const R _phi = arg(ra,ia);
+
+	        rv = _abs*cos(_phi),iv = _abs*sin(_phi);
+        } 
+    };
+
+    // extra
+
+    template<class T> class f_fix {
+    public:
+        /*! \brief Bashes denormals and NANs to zero
+
+	        \param a argument list 
+	        \param v destination vasp (NULL for in-place operation)
+	        \return normalized destination vasp
+        */
+        static V run(T &v,T a) 
+        { 
+	        if(a != a) // NAN
+		        v = 0; 
+	        else {
+		        // denormal bashing (doesn't propagate to the next stage)
+
+		        static const F anti_denormal = 1e-18F;
+		        a += anti_denormal;
+		        a -= anti_denormal;
+		        v = a; 
+	        }
+        } 
+    };
+}
+
+#endif
diff --git a/externals/grill/vasp/source/opparam.cpp b/externals/grill/vasp/source/opparam.cpp
index 6bf42b5e..1c1fd450 100644
--- a/externals/grill/vasp/source/opparam.cpp
+++ b/externals/grill/vasp/source/opparam.cpp
@@ -149,14 +149,35 @@ V OpParam::AI_Rev(I bl)
 		arg[bl].v.idt -= (frames-1)*(arg[bl].v.is = -arg[bl].v.is); 
 }
 
-BL OpParam::AR_In(I bl) const { return arg[bl].argtp == Arg::arg_v && arg[bl].v.rdt && rddt > arg[bl].v.rdt && rddt < arg[bl].v.rdt+frames*arg[bl].v.rs; } 
-BL OpParam::AI_In(I bl) const { return arg[bl].argtp == Arg::arg_v && arg[bl].v.idt && iddt > arg[bl].v.idt && iddt < arg[bl].v.idt+frames*arg[bl].v.is; } 
+BL OpParam::AR_In(I bl) const 
+{ 
+    return arg[bl].argtp == Arg::arg_v && arg[bl].v.rdt && rddt > arg[bl].v.rdt && rddt < arg[bl].v.rdt+frames*arg[bl].v.rs; 
+} 
 
-BL OpParam::AR_Can(I bl) const { return arg[bl].argtp != Arg::arg_v || !arg[bl].v.rdt || arg[bl].v.rdt <= rddt || arg[bl].v.rdt >= rddt+frames*rds; } 
-BL OpParam::AI_Can(I bl) const { return arg[bl].argtp != Arg::arg_v || !arg[bl].v.idt || arg[bl].v.idt <= iddt || arg[bl].v.idt >= iddt+frames*ids; } 
+BL OpParam::AI_In(I bl) const 
+{ 
+    return arg[bl].argtp == Arg::arg_v && arg[bl].v.idt && iddt > arg[bl].v.idt && iddt < arg[bl].v.idt+frames*arg[bl].v.is; 
+} 
 
-BL OpParam::AR_Ovr(I bl) const { return arg[bl].argtp == Arg::arg_v && arg[bl].v.rdt && rddt != arg[bl].v.rdt && rddt < arg[bl].v.rdt+frames*arg[bl].v.rs && arg[bl].v.rdt < rddt+frames*rds; } 
-BL OpParam::AI_Ovr(I bl) const { return arg[bl].argtp == Arg::arg_v && arg[bl].v.idt && iddt != arg[bl].v.idt && iddt < arg[bl].v.idt+frames*arg[bl].v.is && arg[bl].v.idt < iddt+frames*ids; } 
+BL OpParam::AR_Can(I bl) const 
+{ 
+    return arg[bl].argtp != Arg::arg_v || !arg[bl].v.rdt || arg[bl].v.rdt <= rddt || arg[bl].v.rdt >= rddt+frames*rds; 
+} 
+
+BL OpParam::AI_Can(I bl) const 
+{ 
+    return arg[bl].argtp != Arg::arg_v || !arg[bl].v.idt || arg[bl].v.idt <= iddt || arg[bl].v.idt >= iddt+frames*ids; 
+} 
+
+BL OpParam::AR_Ovr(I bl) const 
+{ 
+    return arg[bl].argtp == Arg::arg_v && arg[bl].v.rdt && rddt != arg[bl].v.rdt && rddt < arg[bl].v.rdt+frames*arg[bl].v.rs && arg[bl].v.rdt < rddt+frames*rds; 
+} 
+
+BL OpParam::AI_Ovr(I bl) const 
+{ 
+    return arg[bl].argtp == Arg::arg_v && arg[bl].v.idt && iddt != arg[bl].v.idt && iddt < arg[bl].v.idt+frames*arg[bl].v.is && arg[bl].v.idt < iddt+frames*ids; 
+} 
 
 
 
@@ -170,7 +191,7 @@ BL OpParam::AR_In() const
 BL OpParam::AI_In() const
 {
 	for(I i = 0; i < args; ++i) 
-		if(!AI_In(i)) return true;
+		if(AI_In(i)) return true;
 	return false;
 }
 
diff --git a/externals/grill/vasp/source/ops_arith.cpp b/externals/grill/vasp/source/ops_arith.cpp
index 6faa1aad..5eb11e7d 100644
--- a/externals/grill/vasp/source/ops_arith.cpp
+++ b/externals/grill/vasp/source/ops_arith.cpp
@@ -14,22 +14,6 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 #include "util.h"
 #include <math.h>
 
-template<class T> inline V f_radd(T &v,T a,T b) { v = a+b; }
-template<class T> inline V f_rsub(T &v,T a,T b) { v = a-b; }
-template<class T> inline V f_rsubr(T &v,T a,T b) { v = b-a; }
-template<class T> inline V f_rmul(T &v,T a,T b) { v = a*b; }
-template<class T> inline V f_rdiv(T &v,T a,T b) { v = a/b; }
-template<class T> inline V f_rdivr(T &v,T a,T b) { v = b/a; }
-template<class T> inline V f_rmod(T &v,T a,T b) { v = fmod(a,b); }
-
-BL VecOp::d_add(OpParam &p) { D__rbin(f_radd<S>,p); }
-BL VecOp::d_sub(OpParam &p) { D__rbin(f_rsub<S>,p); }
-BL VecOp::d_subr(OpParam &p) { D__rbin(f_rsubr<S>,p); }
-BL VecOp::d_mul(OpParam &p) { D__rbin(f_rmul<S>,p); }
-BL VecOp::d_div(OpParam &p) { D__rbin(f_rdiv<S>,p); }
-BL VecOp::d_divr(OpParam &p) { D__rbin(f_rdivr<S>,p); }
-BL VecOp::d_mod(OpParam &p) { D__rbin(f_rmod<S>,p); }
-
 
 VASP_BINARY("vasp.+",add,true,VASP_ARG_R(0),"Adds a value, envelope or vasp")
 VASP_BINARY("vasp.-",sub,true,VASP_ARG_R(0),"Subtracts a value, envelope or vasp")
@@ -41,25 +25,11 @@ VASP_BINARY("vasp.%",mod,true,VASP_ARG_R(0),"Calculates the remainder of the div
 
 // -----------------------------------------------------
 
-template<class T> inline V f_rsqr(T &v,T a) { v = a*a; } 
-template<class T> inline V f_rssqr(T &v,T a) { v = a*fabs(a); } 
-
-BL VecOp::d_sqr(OpParam &p) { D__run(f_rsqr<S>,p); }
-BL VecOp::d_ssqr(OpParam &p) { d__run(f_rssqr<S>,p); }
-
 VASP_UNARY("vasp.sqr",sqr,true,"Calculates the square") 
 VASP_UNARY("vasp.ssqr",ssqr,true,"Calculates the square with preservation of the sign") 
 
-
 // -----------------------------------------------------
 
-template<class T> inline V f_rsign(T &v,T a) { v = (a == 0?0:(a < 0?-1.:1.)); }  
-template<class T> inline V f_rabs(T &v,T a) { v = fabs(a); }  
-
-BL VecOp::d_sign(OpParam &p) { D__run(f_rsign<S>,p); }
-BL VecOp::d_abs(OpParam &p) { D__run(f_rabs<S>,p); }
-
-
 VASP_UNARY("vasp.sign",sign,true,"Calculates the sign (signum function)") 
 VASP_UNARY("vasp.abs",abs,true,"Calulates the absolute value") 
 
diff --git a/externals/grill/vasp/source/ops_arith.h b/externals/grill/vasp/source/ops_arith.h
index d4fcd470..7f22b3f6 100644
--- a/externals/grill/vasp/source/ops_arith.h
+++ b/externals/grill/vasp/source/ops_arith.h
@@ -11,34 +11,24 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 #ifndef __VASP_OPS_ARITH_H
 #define __VASP_OPS_ARITH_H
 
-#include "opbase.h"
+#include "opfuns.h"
 
 // Arithmetic math functions
 
 namespace VecOp {
-	BL d_add(OpParam &p); 
-	BL d_sub(OpParam &p); 
-	BL d_subr(OpParam &p); 
-	BL d_mul(OpParam &p); 
-	BL d_div(OpParam &p); 
-	BL d_divr(OpParam &p); 
-	BL d_mod(OpParam &p); 
-
-	BL d_cadd(OpParam &p); 
-	BL d_csub(OpParam &p); 
-	BL d_csubr(OpParam &p); 
-	BL d_cmul(OpParam &p); 
-	BL d_cdiv(OpParam &p); 
-	BL d_cdivr(OpParam &p);
-	
-	BL d_sqr(OpParam &p); 
-	BL d_ssqr(OpParam &p); 
-	BL d_csqr(OpParam &p); 
-	BL d_cpowi(OpParam &p); 
-
-	BL d_sign(OpParam &p); 
-	BL d_abs(OpParam &p); 
-	BL d_cabs(OpParam &p); 
+    inline BL d_add(OpParam &p) { return D__rbin<S,f_add<S> >(p); }
+    inline BL d_sub(OpParam &p) { return D__rbin<S,f_sub<S> >(p); }
+    inline BL d_subr(OpParam &p) { return D__rbin<S,f_subr<S> >(p); }
+    inline BL d_mul(OpParam &p) { return D__rbin<S,f_mul<S> >(p); }
+    inline BL d_div(OpParam &p) { return D__rbin<S,f_div<S> >(p); }
+    inline BL d_divr(OpParam &p) { return D__rbin<S,f_divr<S> >(p); }
+    inline BL d_mod(OpParam &p) { return D__rbin<S,f_mod<S> >(p); }
+
+    inline BL d_sqr(OpParam &p) { return D__run<S,f_sqr<S> >(p); }
+    inline BL d_ssqr(OpParam &p) { return d__run<S,f_ssqr<S> >(p); }
+
+    inline BL d_sign(OpParam &p) { return D__run<S,f_sign<S> >(p); }
+    inline BL d_abs(OpParam &p) { return D__run<S,f_abs<S> >(p); }
 }
 
 namespace VaspOp {
@@ -50,23 +40,11 @@ namespace VaspOp {
 	inline Vasp *m_divr(OpParam &p,CVasp &src,const Argument &arg,CVasp *dst = NULL) { return m_rbin(p,src,arg,dst,VecOp::d_divr); } // reverse div by (one vec or real)
 	inline Vasp *m_mod(OpParam &p,CVasp &src,const Argument &arg,CVasp *dst = NULL) { return m_rbin(p,src,arg,dst,VecOp::d_mod); } // modulo by (one vec or real)
 
-	inline Vasp *m_cadd(OpParam &p,CVasp &src,const Argument &arg,CVasp *dst = NULL) { return m_cbin(p,src,arg,dst,VecOp::d_cadd); }  // complex add (pairs of vecs or complex)
-	inline Vasp *m_csub(OpParam &p,CVasp &src,const Argument &arg,CVasp *dst = NULL) { return m_cbin(p,src,arg,dst,VecOp::d_csub); }  // complex sub (pairs of vecs or complex)
-	inline Vasp *m_csubr(OpParam &p,CVasp &src,const Argument &arg,CVasp *dst = NULL) { return m_cbin(p,src,arg,dst,VecOp::d_csubr); }  // reverse complex sub (pairs of vecs or complex)
-	inline Vasp *m_cmul(OpParam &p,CVasp &src,const Argument &arg,CVasp *dst = NULL) { return m_cbin(p,src,arg,dst,VecOp::d_cmul); }  // complex mul (pairs of vecs or complex)
-	inline Vasp *m_cdiv(OpParam &p,CVasp &src,const Argument &arg,CVasp *dst = NULL) { return m_cbin(p,src,arg,dst,VecOp::d_cdiv); }  // complex div (pairs of vecs or complex)
-	inline Vasp *m_cdivr(OpParam &p,CVasp &src,const Argument &arg,CVasp *dst = NULL) { return m_cbin(p,src,arg,dst,VecOp::d_cdivr); }  // complex reverse div (pairs of vecs or complex)
-
 	inline Vasp *m_sqr(OpParam &p,CVasp &src,CVasp *dst = NULL) { return m_run(p,src,dst,VecOp::d_sqr); }    // unsigned square 
 	inline Vasp *m_ssqr(OpParam &p,CVasp &src,CVasp *dst = NULL) { return m_run(p,src,dst,VecOp::d_ssqr); }   // signed square 
-	inline Vasp *m_csqr(OpParam &p,CVasp &src,CVasp *dst = NULL) { return m_cun(p,src,dst,VecOp::d_csqr); }  // complex square (with each two channels)
-//	inline Vasp *m_csqrt(OpParam &p,CVasp &src,CVasp *dst = NULL) { return m_cun(p,src,dst,VecOp::d_csqrt); }  // complex square root (how about branches?)
-
-	Vasp *m_cpowi(OpParam &p,CVasp &src,const Argument &arg,CVasp *dst = NULL); // complex integer power (with each two channels)
 
 	inline Vasp *m_sign(OpParam &p,CVasp &src,CVasp *dst = NULL) { return m_run(p,src,dst,VecOp::d_sign); }  // sign function 
 	inline Vasp *m_abs(OpParam &p,CVasp &src,CVasp *dst = NULL) { return m_run(p,src,dst,VecOp::d_abs); }  // absolute values
-	inline Vasp *m_cabs(OpParam &p,CVasp &src,CVasp *dst = NULL) { return m_cun(p,src,dst,VecOp::d_cabs); }  // absolute values
 }
 
 #endif
diff --git a/externals/grill/vasp/source/ops_assign.cpp b/externals/grill/vasp/source/ops_assign.cpp
index 05d40375..1d08a373 100644
--- a/externals/grill/vasp/source/ops_assign.cpp
+++ b/externals/grill/vasp/source/ops_assign.cpp
@@ -13,20 +13,6 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 #include "opdefs.h"
 
 
-template<class T> inline V f_rcopy(T &v,T a) { v = a; }
-template<class T> inline V f_ccopy(T &rv,T &iv,T ra,T ia) { rv = ra,iv = ia; }
-
-template<class T> inline V f_rset(T &v,T,T b) { v = b; }
-template<class T> inline V f_cset(T &rv,T &iv,T,T,T rb,T ib) { rv = rb,iv = ib; }
-
-
-BL VecOp::d_copy(OpParam &p) { D__run(f_rcopy<S>,p); }
-BL VecOp::d_ccopy(OpParam &p) { D__cun(f_ccopy<S>,p); }
-
-BL VecOp::d_set(OpParam &p) { D__rbin(f_rset<S>,p); }
-BL VecOp::d_cset(OpParam &p) { D__cbin(f_cset<S>,p); }
-
-
 VASP_BINARY("vasp.set vasp.=",set,false,VASP_ARG_R(0),"Assigns a value, envelope or vasp")
 VASP_BINARY("vasp.cset vasp.c=",cset,false,VASP_ARG_R(0),"Assigns a complex value, real envelope or vasp")
 
diff --git a/externals/grill/vasp/source/ops_assign.h b/externals/grill/vasp/source/ops_assign.h
index e27b3798..a495eb3b 100644
--- a/externals/grill/vasp/source/ops_assign.h
+++ b/externals/grill/vasp/source/ops_assign.h
@@ -11,16 +11,16 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 #ifndef __VASP_OPS_ASSIGN_H
 #define __VASP_OPS_ASSIGN_H
 
-#include "opbase.h"
+#include "opfuns.h"
 
 // Assignment functions
 
 namespace VecOp {
-	BL d_copy(OpParam &p); 
-	BL d_set(OpParam &p); 
+    inline BL d_copy(OpParam &p) { return D__run<S,f_copy<S> >(p); }
+    inline BL d_ccopy(OpParam &p) { return D__cun<S,f_copy<S> >(p); }
 
-	BL d_ccopy(OpParam &p); 
-	BL d_cset(OpParam &p); 
+    inline BL d_set(OpParam &p) { return D__rbin<S,f_set<S> >(p); }
+    inline BL d_cset(OpParam &p) { return D__cbin<S,f_set<S> >(p); }
 }
 
 namespace VaspOp {
diff --git a/externals/grill/vasp/source/ops_carith.cpp b/externals/grill/vasp/source/ops_carith.cpp
index cece7031..5a964cca 100644
--- a/externals/grill/vasp/source/ops_carith.cpp
+++ b/externals/grill/vasp/source/ops_carith.cpp
@@ -15,32 +15,6 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 #include "util.h"
 #include <math.h>
 
-template<class T> inline V f_cadd(T &rv,T &iv,T ra,T ia,T rb,T ib) { rv = ra+rb,iv = ia+ib; }
-template<class T> inline V f_csub(T &rv,T &iv,T ra,T ia,T rb,T ib) { rv = ra-rb,iv = ia-ib; }
-template<class T> inline V f_csubr(T &rv,T &iv,T ra,T ia,T rb,T ib) { rv = rb-ra,iv = ib-ia; }
-template<class T> inline V f_cmul(T &rv,T &iv,T ra,T ia,T rb,T ib) { rv = ra*rb-ia*ib, iv = ra*ib+rb*ia; }
-
-template<class T> inline V f_cdiv(T &rv,T &iv,T ra,T ia,T rb,T ib) 
-{ 
-	register const R den = sqabs(rb,ib);
-	rv = (ra*rb+ia*ib)/den;
-	iv = (ia*rb-ra*ib)/den;
-}
-
-template<class T> inline V f_cdivr(T &rv,T &iv,T ra,T ia,T rb,T ib)
-{ 
-	register const R den = sqabs(ra,ia);
-	rv = (rb*ra+ib*ia)/den;
-	iv = (ib*ra-rb*ia)/den;
-}
-
-BL VecOp::d_cadd(OpParam &p) { D__cbin(f_cadd<S>,p); }
-BL VecOp::d_csub(OpParam &p) { D__cbin(f_csub<S>,p); }
-BL VecOp::d_csubr(OpParam &p) { D__cbin(f_csubr<S>,p); }
-BL VecOp::d_cmul(OpParam &p) { D__cbin(f_cmul<S>,p); }
-BL VecOp::d_cdiv(OpParam &p) { d__cbin(f_cdiv<S>,p); }
-BL VecOp::d_cdivr(OpParam &p) { d__cbin(f_cdivr<S>,p); }
-
 
 VASP_BINARY("vasp.c+",cadd,true,VASP_ARG_R(0),"adds a complex value or vasp")
 VASP_BINARY("vasp.c-",csub,true,VASP_ARG_R(0),"subtracts a complex value or vasp")
@@ -53,24 +27,10 @@ VASP_BINARY("vasp.c!/",cdivr,true,VASP_ARG_R(1),"reverse divides by a complex va
 // -----------------------------------------------------
 
 
-template<class T> inline V f_csqr(T &rv,T &iv,T ra,T ia) { rv = ra*ra-ia*ia; iv = ra*ia*2; }
-
-BL VecOp::d_csqr(OpParam &p) { D__cun(f_csqr<S>,p); }
-
 VASP_UNARY("vasp.csqr",csqr,true,"complex square") 
 
 // -----------------------------------------------------
 
-template<class T> V f_cpowi(T &rv,T &iv,T ra,T ia,OpParam &p) 
-{ 
-	register const I powi = p.ibin.arg;
-	register S rt,it; f_csqr(rt,it,ra,ia);
-	for(I i = 2; i < powi; ++i) f_cmul(rt,it,rt,it,ra,ia);
-	rv = rt,iv = it;
-} 
-
-BL VecOp::d_cpowi(OpParam &p) { d__cop(f_cpowi<S>,p); }
-
 Vasp *VaspOp::m_cpowi(OpParam &p,CVasp &src,const Argument &arg,CVasp *dst) 
 { 
 	Vasp *ret = NULL;
@@ -118,9 +78,5 @@ VASP_ANYOP("vasp.cpowi",cpowi,0,true,VASP_ARG_I(1),"complex integer power")
 
 // -----------------------------------------------------
 
-template<class T> inline V f_cabs(T &rv,T &iv,T ra,T ia) { rv = sqrt(ra*ra+ia*ia),iv = 0; }
-
-BL VecOp::d_cabs(OpParam &p) { D__cun(f_cabs<S>,p); }
-
 VASP_UNARY("vasp.cabs",cabs,true,"set real part to complex absolute value, imaginary part becomes zero") 
 
diff --git a/externals/grill/vasp/source/ops_carith.h b/externals/grill/vasp/source/ops_carith.h
index 3910c0ad..0b049fb4 100644
--- a/externals/grill/vasp/source/ops_carith.h
+++ b/externals/grill/vasp/source/ops_carith.h
@@ -11,22 +11,23 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 #ifndef __VASP_OPS_CARITH_H
 #define __VASP_OPS_CARITH_H
 
-#include "opbase.h"
+#include "opfuns.h"
 
 // Arithmetic math functions
 
 namespace VecOp {
-	BL d_cadd(OpParam &p); 
-	BL d_csub(OpParam &p); 
-	BL d_csubr(OpParam &p); 
-	BL d_cmul(OpParam &p); 
-	BL d_cdiv(OpParam &p); 
-	BL d_cdivr(OpParam &p);
-	
-	BL d_csqr(OpParam &p); 
-	BL d_cpowi(OpParam &p); 
-
-	BL d_cabs(OpParam &p); 
+    inline BL d_cadd(OpParam &p) { return D__cbin<S,f_add<S> >(p); }
+    inline BL d_csub(OpParam &p) { return D__cbin<S,f_sub<S> >(p); }
+    inline BL d_csubr(OpParam &p) { return D__cbin<S,f_subr<S> >(p); }
+    inline BL d_cmul(OpParam &p) { return D__cbin<S,f_mul<S> >(p); }
+    inline BL d_cdiv(OpParam &p) { return d__cbin<S,f_div<S> >(p); }
+    inline BL d_cdivr(OpParam &p) { return d__cbin<S,f_divr<S> >(p); }
+
+    inline BL d_csqr(OpParam &p) { return D__cun<S,f_sqr<S> >(p); }
+
+    inline BL d_cpowi(OpParam &p) { return d__cop<S,f_powi<S> >(p); }
+
+    inline BL d_cabs(OpParam &p) { return D__cun<S,f_abs<S> >(p); }
 }
 
 namespace VaspOp {
diff --git a/externals/grill/vasp/source/ops_cmp.cpp b/externals/grill/vasp/source/ops_cmp.cpp
index fe4bd60a..f2c2ac4d 100644
--- a/externals/grill/vasp/source/ops_cmp.cpp
+++ b/externals/grill/vasp/source/ops_cmp.cpp
@@ -16,28 +16,6 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 
 // --------------------------------------------------------------
 
-template<class T> inline V f_rlwr(T &v,T a,T b) { v = a < b?1:0; }
-template<class T> inline V f_rgtr(T &v,T a,T b) { v = a > b?1:0; }
-template<class T> inline V f_ralwr(T &v,T a,T b) { v = fabs(a) < fabs(b)?1:0; }
-template<class T> inline V f_ragtr(T &v,T a,T b) { v = fabs(a) > fabs(b)?1:0; }
-template<class T> inline V f_rleq(T &v,T a,T b) { v = a <= b?1:0; }
-template<class T> inline V f_rgeq(T &v,T a,T b) { v = a >= b?1:0; }
-template<class T> inline V f_raleq(T &v,T a,T b) { v = fabs(a) <= fabs(b)?1:0; }
-template<class T> inline V f_rageq(T &v,T a,T b) { v = fabs(a) >= fabs(b)?1:0; }
-template<class T> inline V f_requ(T &v,T a,T b) { v = a == b?1:0; }
-template<class T> inline V f_rneq(T &v,T a,T b) { v = a != b?1:0; }
-
-BL VecOp::d_lwr(OpParam &p) { D__rbin(f_rlwr<S>,p); }
-BL VecOp::d_gtr(OpParam &p) { D__rbin(f_rgtr<S>,p); }
-BL VecOp::d_alwr(OpParam &p) { D__rbin(f_ralwr<S>,p); }
-BL VecOp::d_agtr(OpParam &p) { D__rbin(f_ragtr<S>,p); }
-BL VecOp::d_leq(OpParam &p) { D__rbin(f_rleq<S>,p); }
-BL VecOp::d_geq(OpParam &p) { D__rbin(f_rgeq<S>,p); }
-BL VecOp::d_aleq(OpParam &p) { D__rbin(f_raleq<S>,p); }
-BL VecOp::d_ageq(OpParam &p) { D__rbin(f_rageq<S>,p); }
-BL VecOp::d_equ(OpParam &p) { D__rbin(f_requ<S>,p); }
-BL VecOp::d_neq(OpParam &p) { D__rbin(f_rneq<S>,p); }
-
 VASP_BINARY("vasp.<",lwr,true,VASP_ARG_R(0),"set destination to 1 if source < argument, 0 otherwise")
 VASP_BINARY("vasp.>",gtr,true,VASP_ARG_R(0),"set destination to 1 if source > argument, 0 otherwise")
 VASP_BINARY("vasp.a<",alwr,true,VASP_ARG_R(0),"set destination to 1 if abs(source) < abs(argument), 0 otherwise")
@@ -52,27 +30,6 @@ VASP_BINARY("vasp.!=",neq,true,VASP_ARG_R(0),"set destination to 1 if source !=
 
 // --------------------------------------------------------------
 
-template<class T> inline V f_min(T &v,T a,T b) { v = a < b?a:b; }
-template<class T> inline V f_max(T &v,T a,T b) { v = a > b?a:b; }
-
-template<class T> inline V f_rmin(T &rv,T &iv,T ra,T ia,T rb,T ib) 
-{ 
-	if(sqabs(ra,ia) < sqabs(rb,ib))	rv = ra,iv = ia; 
-	else rv = rb,iv = ib; 
-}
-
-template<class T> inline V f_rmax(T &rv,T &iv,T ra,T ia,T rb,T ib) 
-{ 
-	if(sqabs(ra,ia) > sqabs(rb,ib))	rv = ra,iv = ia; 
-	else rv = rb,iv = ib; 
-}
-
-BL VecOp::d_min(OpParam &p) { D__rbin(f_min<S>,p); }
-BL VecOp::d_max(OpParam &p) { D__rbin(f_max<S>,p); }
-BL VecOp::d_rmin(OpParam &p) { d__cbin(f_rmin<S>,p); }
-BL VecOp::d_rmax(OpParam &p) { d__cbin(f_rmax<S>,p); }
-
-
 VASP_BINARY("vasp.min",min,true,VASP_ARG_R(0),"assigns the minimum of the comparison with a value or vasp")
 VASP_BINARY("vasp.max",max,true,VASP_ARG_R(0),"assigns the maximum of the comparison with a value or vasp")
 
@@ -82,14 +39,6 @@ VASP_BINARY("vasp.rmax",rmax,true,VASP_ARG_R(0),"assigns the maximum of the radi
 
 // --------------------------------------------------------------
 
-template<class T> inline V f_minmax(T &rv,T &iv,T ra,T ia) 
-{ 
-	if(ra < ia)	rv = ra,iv = ia; 
-	else rv = ia,iv = ra; 
-} 
-
-BL VecOp::d_minmax(OpParam &p) { d__cun(f_minmax<S>,p); }
-
 VASP_UNARY("vasp, vasp.minmax",minmax,true,"compare two vectors, assign the lower values to the first and the higher to the second one") 
 
 
diff --git a/externals/grill/vasp/source/ops_cmp.h b/externals/grill/vasp/source/ops_cmp.h
index 8fab0e27..73ab132d 100644
--- a/externals/grill/vasp/source/ops_cmp.h
+++ b/externals/grill/vasp/source/ops_cmp.h
@@ -11,42 +11,40 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 #ifndef __VASP_OPS_CMP_H
 #define __VASP_OPS_CMP_H
 
-#include "opbase.h"
+#include "opfuns.h"
 
 // Comparison functions
 
 namespace VecOp {
-	BL d_lwr(OpParam &p); 
-	BL d_gtr(OpParam &p); 
-	BL d_alwr(OpParam &p); 
-	BL d_agtr(OpParam &p); 
-	BL d_leq(OpParam &p); 
-	BL d_geq(OpParam &p); 
-	BL d_aleq(OpParam &p); 
-	BL d_ageq(OpParam &p); 
-	BL d_equ(OpParam &p); 
-	BL d_neq(OpParam &p); 
-
-	BL d_min(OpParam &p); 
-	BL d_max(OpParam &p); 
-
-	BL d_rmin(OpParam &p); 
-	BL d_rmax(OpParam &p); 
-
-	BL d_minmax(OpParam &p); 
-
-	BL d_minq(OpParam &p); 
-	BL d_maxq(OpParam &p); 
-	BL d_aminq(OpParam &p); 
-	BL d_amaxq(OpParam &p); 
-
-	BL d_rminq(OpParam &p); 
-	BL d_rmaxq(OpParam &p); 
-
-	BL d_gate(OpParam &p); 
-	BL d_igate(OpParam &p); 
-	BL d_rgate(OpParam &p); 
-	BL d_rigate(OpParam &p); 
+    inline BL d_lwr(OpParam &p) { return D__rbin<S,f_lwr<S> >(p); }
+    inline BL d_gtr(OpParam &p) { return D__rbin<S,f_gtr<S> >(p); }
+    inline BL d_alwr(OpParam &p) { return D__rbin<S,f_alwr<S> >(p); }
+    inline BL d_agtr(OpParam &p) { return D__rbin<S,f_agtr<S> >(p); }
+    inline BL d_leq(OpParam &p) { return D__rbin<S,f_leq<S> >(p); }
+    inline BL d_geq(OpParam &p) { return D__rbin<S,f_geq<S> >(p); }
+    inline BL d_aleq(OpParam &p) { return D__rbin<S,f_aleq<S> >(p); }
+    inline BL d_ageq(OpParam &p) { return D__rbin<S,f_ageq<S> >(p); }
+    inline BL d_equ(OpParam &p) { return D__rbin<S,f_equ<S> >(p); }
+    inline BL d_neq(OpParam &p) { return D__rbin<S,f_neq<S> >(p); }
+
+    inline BL d_min(OpParam &p) { return D__rbin<S,f_min<S> >(p); }
+    inline BL d_max(OpParam &p) { return D__rbin<S,f_max<S> >(p); }
+    inline BL d_rmin(OpParam &p) { return d__cbin<S,f_min<S> >(p); }
+    inline BL d_rmax(OpParam &p) { return d__cbin<S,f_max<S> >(p); }
+
+    inline BL d_minmax(OpParam &p) { return d__cun<S,f_minmax<S> >(p); }
+
+    inline BL d_minq(OpParam &p) { return D__rop<S,f_minq<S> >(p); }
+    inline BL d_maxq(OpParam &p) { return D__rop<S,f_maxq<S> >(p); }
+    inline BL d_rminq(OpParam &p) { return d__cop<S,f_minq<S> >(p); }
+    inline BL d_rmaxq(OpParam &p) { return d__cop<S,f_maxq<S> >(p); }
+    inline BL d_aminq(OpParam &p) { return F__rop<S,f_aminq<S> >(p); }
+    inline BL d_amaxq(OpParam &p) { return F__rop<S,f_amaxq<S> >(p); }
+
+    inline BL d_gate(OpParam &p) { return D__rbin<S,f_gate<S> >(p); }
+    inline BL d_igate(OpParam &p) { return d__rbin<S,f_igate<S> >(p); }
+    inline BL d_rgate(OpParam &p) { return d__cbin<S,f_gate<S> >(p); }
+    inline BL d_rigate(OpParam &p) { return d__cbin<S,f_igate<S> >(p); }
 }
 
 namespace VaspOp {
diff --git a/externals/grill/vasp/source/ops_cplx.cpp b/externals/grill/vasp/source/ops_cplx.cpp
index b2f01788..013571ef 100644
--- a/externals/grill/vasp/source/ops_cplx.cpp
+++ b/externals/grill/vasp/source/ops_cplx.cpp
@@ -16,13 +16,6 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 
 // -----------------------------------------------------
 
-template<class T> V f_polar(T &rv,T &iv,T ra,T ia) { rv = sqrt(sqabs(ra,ia)),iv = arg(ra,ia); }
-template<class T> V f_rect(T &rv,T &iv,T ra,T ia) { rv = ra*cos(ia),iv = ra*sin(ia); }
-
-BL VecOp::d_polar(OpParam &p) { d__cun(f_polar<S>,p); }
-BL VecOp::d_rect(OpParam &p) { d__cun(f_rect<S>,p); }
-
-
 VASP_UNARY("vasp.polar",polar,true,"convert complex vector pair from rectangular to polar coordinates") 
 VASP_UNARY("vasp.rect",rect,true,"convert complex vector pair from polar to rectangular coordinates") 
 
@@ -30,16 +23,6 @@ VASP_UNARY("vasp.rect",rect,true,"convert complex vector pair from polar to rect
 // -----------------------------------------------------
 
 
-template<class T> V f_radd(T &rv,T &iv,T ra,T ia,T rb,T) 
-{ 
-	register const R _abs = sqrt(sqabs(ra,ia))+rb;
-	register const R _phi = arg(ra,ia);
-
-	rv = _abs*cos(_phi),iv = _abs*sin(_phi);
-} 
-
-BL VecOp::d_radd(OpParam &p) { d__cbin(f_radd<S>,p); }
-
 Vasp *VaspOp::m_radd(OpParam &p,CVasp &src,const Argument &arg,CVasp *dst) 
 { 
 	Vasp *ret = NULL;
@@ -65,22 +48,9 @@ VASP_ANYOP("vasp.r+",radd,0,true,VASP_ARG_R(0),"add offset to complex radius (of
 
 // -----------------------------------------------------
 
-template<class T> V f_cnorm(T &rv,T &iv,T ra,T ia) 
-{ 
-	register T f = sqabs(ra,ia);
-	if(f) { f = 1./sqrt(f); rv = ra*f,iv = ia*f; }
-	else rv = iv = 0;
-}
-
-BL VecOp::d_cnorm(OpParam &p) { d__cun(f_cnorm<S>,p); }
-
 VASP_UNARY("vasp.cnorm",cnorm,true,"normalize complex radius to 1 (but preserve angle)")
 
 // -----------------------------------------------------
 
-template<class T> inline V f_cconj(T &,T &iv,T,T ia) { iv = -ia; }
-
-BL VecOp::d_cconj(OpParam &p) { D__cun(f_cconj<S>,p); }
-
 VASP_UNARY("vasp.cconj",cconj,true,"complex conjugate: multiply imaginary part with -1")  // should be replaced by an abstraction
 
diff --git a/externals/grill/vasp/source/ops_cplx.h b/externals/grill/vasp/source/ops_cplx.h
index b3aaa4c6..c77de675 100644
--- a/externals/grill/vasp/source/ops_cplx.h
+++ b/externals/grill/vasp/source/ops_cplx.h
@@ -11,20 +11,19 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 #ifndef __VASP_OPS_CPLX_H
 #define __VASP_OPS_CPLX_H
 
-#include "opbase.h"
+#include "opfuns.h"
 
 // Complex functions
 
 namespace VecOp {
-	BL d_polar(OpParam &p); 
-	BL d_rect(OpParam &p); 
+    inline BL d_polar(OpParam &p) { return d__cun<S,f_polar<S> >(p); }
+    inline BL d_rect(OpParam &p) { return d__cun<S,f_rect<S> >(p); }
 
-	BL d_radd(OpParam &p); 
+    inline BL d_radd(OpParam &p) { return d__cbin<S,f_radd<S> >(p); }
 
-	BL d_cnorm(OpParam &p); 
+    inline BL d_cnorm(OpParam &p) { return d__cun<S,f_norm<S> >(p); }
 
-//	BL d_cswap(OpParam &p); 
-	BL d_cconj(OpParam &p); 
+    inline BL d_cconj(OpParam &p) { return D__cun<S,f_conj<S> >(p); }
 }
 
 namespace VaspOp {
diff --git a/externals/grill/vasp/source/ops_flt.cpp b/externals/grill/vasp/source/ops_flt.cpp
index fae66f56..f92f585e 100644
--- a/externals/grill/vasp/source/ops_flt.cpp
+++ b/externals/grill/vasp/source/ops_flt.cpp
@@ -209,27 +209,4 @@ VASP_ANYOP("vasp.int",int,0,true,VASP_ARG_I(1),"Integration")
 VASP_ANYOP("vasp.dif",dif,0,true,VASP_ARG_I(1),"Differentiation") 
 
 
-
-/*! \brief Bashes denormals and NANs to zero
-
-	\param arg argument list 
-	\param dst destination vasp (NULL for in-place operation)
-	\return normalized destination vasp
-*/
-template<class T> inline V f_fix(T &v,T a) 
-{ 
-	if(a != a) // NAN
-		v = 0; 
-	else {
-		// denormal bashing (doesn't propagate to the next stage)
-
-		static const F anti_denormal = 1e-18F;
-		a += anti_denormal;
-		a -= anti_denormal;
-		v = a; 
-	}
-} 
-
-BL VecOp::d_fix(OpParam &p) { D__run(f_fix<S>,p); }
-
 VASP_UNARY("vasp.fix",fix,true,"Bashes denormals/NANs to zero") 
diff --git a/externals/grill/vasp/source/ops_flt.h b/externals/grill/vasp/source/ops_flt.h
index 21761d90..ea7fe5ce 100644
--- a/externals/grill/vasp/source/ops_flt.h
+++ b/externals/grill/vasp/source/ops_flt.h
@@ -11,7 +11,7 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 #ifndef __VASP_OPS_FLT_H
 #define __VASP_OPS_FLT_H
 
-#include "opbase.h"
+#include "opfuns.h"
 
 // Filtering functions
 
@@ -22,7 +22,7 @@ namespace VecOp {
 	BL d_int(OpParam &p);
 	BL d_dif(OpParam &p); 
 
-	BL d_fix(OpParam &p); 
+    inline BL d_fix(OpParam &p) { return D__run<S,f_fix<S> >(p); }
 }
 
 namespace VaspOp {
diff --git a/externals/grill/vasp/source/ops_gate.cpp b/externals/grill/vasp/source/ops_gate.cpp
index 3dd15fbb..919ee07f 100644
--- a/externals/grill/vasp/source/ops_gate.cpp
+++ b/externals/grill/vasp/source/ops_gate.cpp
@@ -17,32 +17,6 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 // --------------------------------------------------------------
 
 
-template<class T> V f_gate(T &rv,T ra,T rb) { rv = fabs(ra) >= rb?ra:0; } 
-template<class T> V f_igate(T &rv,T ra,T rb) { rv = fabs(ra) <= rb?ra:0; } 
-
-template<class T> V f_rgate(T &rv,T &iv,T ra,T ia,T rb,T) 
-{ 
-	register const T _abs = sqabs(ra,ia);
-
-	if(_abs >= rb*rb) rv = ra,iv = ia;
-	else rv = iv = 0;
-} 
-
-template<class T> V f_rigate(T &rv,T &iv,T ra,T ia,T rb,T) 
-{ 
-	register const T _abs = sqabs(ra,ia);
-
-	if(_abs <= rb*rb) rv = ra,iv = ia;
-	else rv = iv = 0;
-} 
-
-BL VecOp::d_gate(OpParam &p) { D__rbin(f_gate<S>,p); }
-BL VecOp::d_igate(OpParam &p) { d__rbin(f_igate<S>,p); }
-BL VecOp::d_rgate(OpParam &p) { d__cbin(f_rgate<S>,p); }
-BL VecOp::d_rigate(OpParam &p) { d__cbin(f_rigate<S>,p); }
-
-
-
 Vasp *VaspOp::m_gate(OpParam &p,CVasp &src,const Argument &arg,CVasp *dst) 
 { 
 	Vasp *ret = NULL;
diff --git a/externals/grill/vasp/source/ops_qminmax.cpp b/externals/grill/vasp/source/ops_qminmax.cpp
index 6fe8ee26..f7e4144e 100644
--- a/externals/grill/vasp/source/ops_qminmax.cpp
+++ b/externals/grill/vasp/source/ops_qminmax.cpp
@@ -16,47 +16,6 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 
 // --------------------------------------------------------------
 
-template<class T> V f_maxq(T &,T ra,OpParam &p) 
-{ 
-	if(ra > p.norm.minmax) p.norm.minmax = ra; 
-} 
-
-template<class T> V f_minq(T &,T ra,OpParam &p) 
-{ 
-	if(ra < p.norm.minmax) p.norm.minmax = ra; 
-} 
-
-template<class T> V f_amaxq(T &,T ra,OpParam &p) 
-{ 
-	register T s = fabs(ra); 
-	if(s > p.norm.minmax) p.norm.minmax = s; 
-} 
-
-template<class T> V f_aminq(T &,T ra,OpParam &p) 
-{ 
-	register T s = fabs(ra); 
-	if(s < p.norm.minmax) p.norm.minmax = s; 
-} 
-
-template<class T> V f_rmaxq(T &,T &,T ra,T ia,OpParam &p) 
-{ 
-	register T s = sqabs(ra,ia); 
-	if(s > p.norm.minmax) p.norm.minmax = s; 
-} 
-
-template<class T> V f_rminq(T &,T &,T ra,T ia,OpParam &p) 
-{ 
-	register T s = sqabs(ra,ia); 
-	if(s < p.norm.minmax) p.norm.minmax = s; 
-} 
-
-BL VecOp::d_minq(OpParam &p) { D__rop(f_minq<S>,p); }
-BL VecOp::d_maxq(OpParam &p) { D__rop(f_maxq<S>,p); }
-BL VecOp::d_aminq(OpParam &p) { D__rop(f_aminq<S>,p); }
-BL VecOp::d_amaxq(OpParam &p) { D__rop(f_amaxq<S>,p); }
-BL VecOp::d_rminq(OpParam &p) { d__cop(f_rminq<S>,p); }
-BL VecOp::d_rmaxq(OpParam &p) { d__cop(f_rmaxq<S>,p); }
-
 // --------------------------------------------------------------
 
 
diff --git a/externals/grill/vasp/source/ops_trnsc.cpp b/externals/grill/vasp/source/ops_trnsc.cpp
index 88089fbb..8bbda8e9 100644
--- a/externals/grill/vasp/source/ops_trnsc.cpp
+++ b/externals/grill/vasp/source/ops_trnsc.cpp
@@ -16,23 +16,6 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 
 // --------------------------------------------------------------
 
-template<class T> V f_rpow(T &v,T a,T b) { v = pow(fabs(a),b)*sgn(a); } 
-
-BL VecOp::d_pow(OpParam &p) { d__rbin(f_rpow<S>,p); }
-
-template<class T> V f_crpow(T &rv,T &iv,T ra,T ia,T rb,T) 
-{ 
-	register const R _abs = sqrt(sqabs(ra,ia));
-	if(_abs) {
-		register const R _p = pow(_abs,rb)/_abs;
-		rv = _p*ra,iv = _p*ia;
-	}
-	else
-		rv = iv = 0;
-} 
-
-BL VecOp::d_rpow(OpParam &p) { d__cbin(f_crpow<S>,p); }
-
 Vasp *VaspOp::m_rpow(OpParam &p,CVasp &src,const Argument &arg,CVasp *dst) 
 { 
 	Vasp *ret = NULL;
@@ -58,26 +41,11 @@ VASP_ANYOP("vasp.rpow",rpow,0,true,VASP_ARG_R(1),"Power function acting on compl
 
 // --------------------------------------------------------------
 
-template<class T> V f_rsqrt(T &v,T a) { v = sqrt(fabs(a)); } 
-template<class T> V f_rssqrt(T &v,T a) { v = sqrt(fabs(a))*sgn(a); } 
-
-BL VecOp::d_sqrt(OpParam &p) { d__run(f_rsqrt<S>,p); }
-BL VecOp::d_ssqrt(OpParam &p) { d__run(f_rssqrt<S>,p); }
-
-
 VASP_UNARY("vasp.sqrt",sqrt,true,"Square root") 
 VASP_UNARY("vasp.ssqrt",ssqrt,true,"Square root preserving the sign") 
 
 // --------------------------------------------------------------
 
-
-template<class T> V f_rexp(T &v,T a) { v = exp(a); } 
-template<class T> V f_rlog(T &v,T a) { v = log(a); }  // \todo detect NANs
-
-BL VecOp::d_exp(OpParam &p) { d__run(f_rexp<S>,p); }
-BL VecOp::d_log(OpParam &p) { d__run(f_rlog<S>,p); }
-
-
 VASP_UNARY("vasp.exp",exp,true,"Exponential function") 
 VASP_UNARY("vasp.log",log,true,"Natural logarithm") 
 
diff --git a/externals/grill/vasp/source/ops_trnsc.h b/externals/grill/vasp/source/ops_trnsc.h
index ee919ea9..73617d9a 100644
--- a/externals/grill/vasp/source/ops_trnsc.h
+++ b/externals/grill/vasp/source/ops_trnsc.h
@@ -11,21 +11,19 @@ WARRANTIES, see the file, "license.txt," in this distribution.
 #ifndef __VASP_OPS_TRNSC_H
 #define __VASP_OPS_TRNSC_H
 
-#include "opbase.h"
+#include "opfuns.h"
 
 // Transcendent math functions
 
 namespace VecOp {
-	BL d_pow(OpParam &p); 
+    inline BL d_pow(OpParam &p) { return d__rbin<S,f_pow<S> >(p); }
+    inline BL d_rpow(OpParam &p) { return d__cbin<S,f_pow<S> >(p); }
 
-	BL d_rpow(OpParam &p); 
+    inline BL d_sqrt(OpParam &p) { return d__run<S,f_sqrt<S> >(p); }
+    inline BL d_ssqrt(OpParam &p) { return d__run<S,f_ssqrt<S> >(p); }
 
-	BL d_sqrt(OpParam &p); 
-	BL d_ssqrt(OpParam &p); 
-
-
-	BL d_exp(OpParam &p); 
-	BL d_log(OpParam &p); 
+    inline BL d_exp(OpParam &p) { return d__run<S,f_exp<S> >(p); }
+    inline BL d_log(OpParam &p) { return d__run<S,f_log<S> >(p); }
 }
 
 namespace VaspOp {
diff --git a/externals/grill/vasp/vasp.dsp b/externals/grill/vasp/vasp.dsp
index d3077050..9fbbd71a 100644
--- a/externals/grill/vasp/vasp.dsp
+++ b/externals/grill/vasp/vasp.dsp
@@ -377,6 +377,10 @@ SOURCE=source\opdefs.h
 # End Source File
 # Begin Source File
 
+SOURCE=source\opfuns.h
+# End Source File
+# Begin Source File
+
 SOURCE=source\oploop.h
 # End Source File
 # Begin Source File