diff --git a/mex/sources/kalman/qt/FortranSubroutines.zip b/mex/sources/kalman/qt/FortranSubroutines.zip
new file mode 100644
index 000000000..5a69ed6b9
Binary files /dev/null and b/mex/sources/kalman/qt/FortranSubroutines.zip differ
diff --git a/mex/sources/kalman/qt/cc/qt.h b/mex/sources/kalman/qt/cc/qt.h
new file mode 100644
index 000000000..ed4637bac
--- /dev/null
+++ b/mex/sources/kalman/qt/cc/qt.h
@@ -0,0 +1,32 @@
+
+#ifndef QT_H
+#define QT_H
+
+#define C_CHAR const char*
+#define BLINT  int*
+#define C_BLINT const int*
+#define C_BLDOU const double*
+#define BLDOU double*
+
+
+extern "C"{
+void ldm_(BLDOU, C_BLDOU, C_BLDOU, C_BLINT);
+void ldsld_(BLDOU, C_BLDOU, C_BLDOU, C_BLINT);
+void ldv_(BLDOU,C_BLDOU, C_BLDOU, C_BLINT);
+void mtt_(BLDOU,C_BLDOU, C_BLDOU, C_BLINT);
+void qt2ld_(BLDOU,C_BLDOU, C_BLINT);
+void qt2t_(BLDOU, C_BLDOU, C_BLINT);
+void s2d_(BLDOU,C_BLDOU, C_BLINT);
+void s2u_(BLDOU,C_BLDOU, C_BLINT);
+void td_(BLDOU,C_BLDOU, C_BLDOU, C_BLINT);
+void tm_(BLDOU,C_BLDOU, C_BLDOU, C_BLINT);
+void tstt_(BLDOU,C_BLDOU, C_BLDOU, C_BLINT);
+void tt_(BLDOU,C_BLDOU, C_BLDOU, C_BLINT);
+void tu_(BLDOU,C_BLDOU, C_BLDOU, C_BLINT);
+void tut_(BLDOU,C_BLDOU, C_BLDOU, C_BLINT);
+void tv_(BLDOU,C_BLDOU, C_BLDOU, C_BLINT);
+void qtv_(BLDOU,C_BLDOU, C_BLDOU, C_BLINT);
+void qtsqtt_(BLDOU,C_BLDOU, C_BLDOU, C_BLINT);
+
+};
+#endif
diff --git a/mex/sources/kalman/qt/f90/LdM.f90 b/mex/sources/kalman/qt/f90/LdM.f90
new file mode 100644
index 000000000..8f9f1a765
--- /dev/null
+++ b/mex/sources/kalman/qt/f90/LdM.f90
@@ -0,0 +1,24 @@
+subroutine LdM(X,Ld,M,n)
+!     COMPUTATIONAL SUBROUTINE: LdM=Ld*M; Ld lower diagonal matrix; M arbitrary
+
+implicit none
+
+integer(4), intent(in) :: n
+real(8), intent(in) :: Ld(n,n)
+real(8), intent(in) :: M(n,n)
+real(8), intent(out) :: X(n,n)
+integer(4) :: i,j
+
+
+X=0.0*M
+
+do i=2,n
+if ((Ld(i,i-1)/=0.0)) then
+do j=1,n
+X(i,j)=Ld(i,i-1)*M(i-1,j)
+end do
+end if
+end do
+
+return
+end subroutine LdM
\ No newline at end of file
diff --git a/mex/sources/kalman/qt/f90/LdSLd.f90 b/mex/sources/kalman/qt/f90/LdSLd.f90
new file mode 100644
index 000000000..1fd861143
--- /dev/null
+++ b/mex/sources/kalman/qt/f90/LdSLd.f90
@@ -0,0 +1,34 @@
+subroutine LdSLd(X,Ld,S,n)
+!     COMPUTATIONAL SUBROUTINE: LdM=Ld*S*Ld; Ld lower diagonal matrix; S symmetric
+
+implicit none
+
+integer(4), intent(in) :: n
+real(8), intent(in) :: Ld(n,n)
+real(8), intent(in) :: S(n,n)
+real(8), intent(out) :: X(n,n)
+integer(4), dimension(1,n) :: vk
+
+
+integer(4) :: i,j, jj, h
+
+
+h=0.0
+do  i=2,n
+if ((Ld(i,i-1)/= 0.0)) then
+h=h+1.0
+vk(1,h)=i
+end if
+end do
+if (h==0.0) then
+X=0*S
+return
+end if
+do  j=1,h
+do  jj=1,h
+X(vk(1,j),vk(1,jj))=Ld(vk(1,j),vk(1,j)-1)*S(vk(1,j)-1,vk(1,jj)-1)*Ld(vk(1,jj),vk(1,jj)-1)
+end do
+end do
+
+return
+end subroutine LdSLd
diff --git a/mex/sources/kalman/qt/f90/LdV.f90 b/mex/sources/kalman/qt/f90/LdV.f90
new file mode 100644
index 000000000..66671a749
--- /dev/null
+++ b/mex/sources/kalman/qt/f90/LdV.f90
@@ -0,0 +1,17 @@
+subroutine LdV(X,Ld,v,n)
+implicit none
+
+integer(4), intent(in) :: n
+real(8), intent(in) :: ld(n,n)
+real(8), intent(in) :: v(n,1)
+real(8), intent(out) :: X(n,1)
+integer(4) :: i
+
+X=0*v
+
+do i=2,n
+X(i,1)=Ld(i,i-1)*v(i-1,1)
+end do
+
+return
+end subroutine LdV
diff --git a/mex/sources/kalman/qt/f90/MTt.f90 b/mex/sources/kalman/qt/f90/MTt.f90
new file mode 100644
index 000000000..54daa7d1a
--- /dev/null
+++ b/mex/sources/kalman/qt/f90/MTt.f90
@@ -0,0 +1,31 @@
+subroutine MTt(X,M,Tt,n)
+!     COMPUTATIONAL SUBROUTINE TUt=T*Ut T upper triangular; U striclty lower triangular
+
+implicit none
+
+integer(4), intent(in) :: n
+real(8), intent(in) :: M(n,n)
+real(8), intent(in) :: Tt(n,n)
+real(8), intent(out) :: X(n,n)
+
+integer(4) :: i,j,k
+
+real(8) :: stemp
+
+
+
+X=0.0*M
+
+do i=1,n
+do j=1,n
+stemp = 0.0
+do  k = j,n
+stemp = stemp + M(i,k+i) * Tt(k+i,j)
+end do
+X(i,j)=stemp
+end do
+
+end do
+
+return
+end subroutine MTt
\ No newline at end of file
diff --git a/mex/sources/kalman/qt/f90/QT2Ld.f90 b/mex/sources/kalman/qt/f90/QT2Ld.f90
new file mode 100644
index 000000000..f542fecad
--- /dev/null
+++ b/mex/sources/kalman/qt/f90/QT2Ld.f90
@@ -0,0 +1,26 @@
+subroutine QT2Ld(X,QT,n)
+!     COMPUTATIONAL SUBROUTINE: extracts lower diagonal from Quasi triangular matrix
+
+!     COMPUTATIONAL SUBROUTINE: extracts lower diagonal from Quasi triangular matrix
+implicit none
+
+integer(4), intent(in) :: n
+real(8), intent(in) :: QT(n,n)
+real(8), intent(out) :: X(n,n)
+integer(4) :: i
+
+
+
+
+
+X=0.0*QT
+
+
+do i=2,n
+if (QT(i,i-1)/=0.0) then
+X(i,i-1)=QT(i,i-1)
+end if
+end do
+
+return
+end subroutine QT2Ld
diff --git a/mex/sources/kalman/qt/f90/QT2T.f90 b/mex/sources/kalman/qt/f90/QT2T.f90
new file mode 100644
index 000000000..0eef455ec
--- /dev/null
+++ b/mex/sources/kalman/qt/f90/QT2T.f90
@@ -0,0 +1,20 @@
+subroutine QT2T(X,QT,n)
+!     COMPUTATIONAL SUBROUTINE: extracts upper triangular from Quasi triangular matrix
+
+implicit none
+
+integer(4), intent(in) :: n
+real(8), intent(in) :: QT(n,n)
+real(8), intent(out) :: X(n,n)
+integer(4) :: i,j
+
+
+X=0.0*QT
+do i=1,n
+do j=i,n
+X(i,j)=QT(i,j)
+end do
+end do
+
+return
+end subroutine QT2T
diff --git a/mex/sources/kalman/qt/f90/QTSQTt.f90 b/mex/sources/kalman/qt/f90/QTSQTt.f90
new file mode 100644
index 000000000..283e55585
--- /dev/null
+++ b/mex/sources/kalman/qt/f90/QTSQTt.f90
@@ -0,0 +1,38 @@
+subroutine QTSQTt(X,QT,S,n)
+!     COMPUTATIONAL SUBROUTINE TUt=T*Ut T upper triangular; U striclty lower triangular
+
+implicit none
+
+integer(4), intent(in) :: n
+real(8), intent(in) :: QT(n,n)
+real(8), intent(in) :: S(n,n)
+
+real(8), intent(out) :: X(n,n)
+integer(4) :: i,j,k,h
+
+real(8) :: stemp
+
+
+X=0.0*S
+
+do i=1,n
+do j=1,n
+stemp = 0.0
+if (i > 1 .AND. (QT(i,i-1)/= 0.0)) then
+stemp = QT(i,i-1)*S(i-1,1)*QT(i,i-1)
+end if 
+
+do  h = i,n
+do  k = j,n
+stemp = stemp + QT(i,h) * S(h,k) * QT(j,k)
+end do
+end do
+X(i,j)=stemp
+X(j,i)=stemp
+
+end do
+
+end do
+
+return
+end subroutine QTSQTt
diff --git a/mex/sources/kalman/qt/f90/QTV.f90 b/mex/sources/kalman/qt/f90/QTV.f90
new file mode 100644
index 000000000..fbe2cc5e4
--- /dev/null
+++ b/mex/sources/kalman/qt/f90/QTV.f90
@@ -0,0 +1,33 @@
+subroutine QTV(X,QT,V,n)
+!     COMPUTATIONAL SUBROUTINE: X=QT*V QT upper quasi-triangular; V vector
+
+implicit none
+
+integer(4), intent(in) :: n
+real(8), intent(in) :: QT(n,n)
+real(8), intent(in) :: V(n,1)
+
+real(8), intent(out) :: X(n,1) 
+
+integer(4) :: i,k
+
+real(8) :: stemp
+
+
+
+do i=1,n
+stemp = 0.0
+if (i > 1 .AND. (QT(i,i-1)/= 0.0)) then
+stemp = QT(i,i-1)*v(i-1,1)
+end if 
+
+do  k = 0,n-i
+stemp = stemp + QT(i,k+i) * V(k+i,1)
+end do
+X(i,1)=stemp
+end do
+
+return
+end subroutine QTV
+
+
diff --git a/mex/sources/kalman/qt/f90/S2D.f90 b/mex/sources/kalman/qt/f90/S2D.f90
new file mode 100644
index 000000000..432a73016
--- /dev/null
+++ b/mex/sources/kalman/qt/f90/S2D.f90
@@ -0,0 +1,19 @@
+subroutine S2D(X,S,n)
+!     COMPUTATIONAL SUBROUTINE TUt=T*Ut T upper triangular; U striclty upper triangular
+
+implicit none
+
+integer(4), intent(in) :: n
+real(8), intent(in) :: S(n,n)
+real(8), intent(out) :: X(n,n)
+
+integer(4) :: i
+
+
+
+X=0.0*S
+do i=1,n
+X(i,i)=sqrt(S(i,i))
+end do
+return
+end subroutine S2D
diff --git a/mex/sources/kalman/qt/f90/S2U.f90 b/mex/sources/kalman/qt/f90/S2U.f90
new file mode 100644
index 000000000..d7ce4dece
--- /dev/null
+++ b/mex/sources/kalman/qt/f90/S2U.f90
@@ -0,0 +1,23 @@
+subroutine S2U(X,S,n)
+!     COMPUTATIONAL SUBROUTINE: S2U extracts striclty upper triangular from symmetric
+
+implicit none
+
+integer(4), intent(in) :: n
+real(8), intent(in) :: S(n,n)
+real(8), intent(out) :: X(n,n)
+integer(4) :: i,j
+
+
+
+
+
+X=0.0*S
+do i=1,n
+do j=i+1,n
+X(i,j)=S(i,j)
+end do
+end do
+
+return
+end subroutine S2U
diff --git a/mex/sources/kalman/qt/f90/TD.f90 b/mex/sources/kalman/qt/f90/TD.f90
new file mode 100644
index 000000000..59dd31f31
--- /dev/null
+++ b/mex/sources/kalman/qt/f90/TD.f90
@@ -0,0 +1,22 @@
+subroutine TD(X,T,D,n)
+!     COMPUTATIONAL SUBROUTINE TUt=T*D T upper triangular; D diagonal
+
+implicit none
+
+integer(4), intent(in) :: n
+real(8), intent(in) :: T(n,n)
+real(8), intent(in) :: D(n,n)
+
+real(8), intent(out) :: X(n,n)
+integer(4) :: i,j
+
+
+X=0.0*T
+do i=1,n
+do j=i,n
+X(i,j)=T(i,j)*D(j,j)
+end do
+end do
+
+return
+end subroutine TD
diff --git a/mex/sources/kalman/qt/f90/TM.f90 b/mex/sources/kalman/qt/f90/TM.f90
new file mode 100644
index 000000000..4e0f2f1fb
--- /dev/null
+++ b/mex/sources/kalman/qt/f90/TM.f90
@@ -0,0 +1,29 @@
+subroutine TM(X,T,M,n)
+!     COMPUTATIONAL SUBROUTINE: TM=T*M T upper triangular; M arbitrary
+
+
+implicit none
+
+integer(4), intent(in) :: n
+real(8), intent(in) :: T(n,n)
+real(8), intent(in) :: M(n,n)
+
+real(8), intent(out) :: X(n,n)
+
+integer(4) :: i,j,k
+
+real(8) :: stemp
+
+
+do i=1,n
+do j=1,n
+stemp = 0.0
+do  k = 0,n-i
+stemp = stemp + T(i,k+i) * M(k+i,j)
+end do
+X(i,j)=stemp
+end do
+end do
+
+return
+end subroutine TM
diff --git a/mex/sources/kalman/qt/f90/TSTt.f90 b/mex/sources/kalman/qt/f90/TSTt.f90
new file mode 100644
index 000000000..d44ffd76d
--- /dev/null
+++ b/mex/sources/kalman/qt/f90/TSTt.f90
@@ -0,0 +1,34 @@
+subroutine TSTt(X,T,S,n)
+!     COMPUTATIONAL SUBROUTINE TUt=T*Ut T upper triangular; U striclty lower triangular
+
+implicit none
+
+integer(4), intent(in) :: n
+real(8), intent(in) :: T(n,n)
+real(8), intent(in) :: S(n,n)
+
+real(8), intent(out) :: X(n,n)
+integer(4) :: i,j,k,h
+
+real(8) :: stemp
+
+
+X=0.0*S
+
+do i=1,n
+do j=1,n
+stemp = 0.0
+do  h = i,n
+do  k = j,n
+stemp = stemp + T(i,h) * S(h,k) * T(j,k)
+end do
+end do
+X(i,j)=stemp
+X(j,i)=stemp
+
+end do
+
+end do
+
+return
+end subroutine TSTt
diff --git a/mex/sources/kalman/qt/f90/TT.f90 b/mex/sources/kalman/qt/f90/TT.f90
new file mode 100644
index 000000000..5024942e0
--- /dev/null
+++ b/mex/sources/kalman/qt/f90/TT.f90
@@ -0,0 +1,29 @@
+subroutine TT(X, T1,T2,n)
+!     COMPUTATIONAL SUBROUTINE TUt=T*Ut T upper triangular; U striclty lower triangular
+
+implicit none
+
+integer(4), intent(in) :: n
+real(8), intent(in) :: T1(n,n)
+real(8), intent(in) :: T2(n,n)
+
+real(8), intent(out) :: X(n,n)
+integer(4) :: i,j,k,h
+
+real(8) :: stemp
+
+X=0.0*T1
+
+do i=1,n
+do j=1,n
+stemp = 0.0
+do  k = 0,n-i
+stemp = stemp + T1(i,k+i) * T2(k+i,j)
+end do
+X(i,j)=stemp
+end do
+
+end do
+
+return
+end subroutine TT
\ No newline at end of file
diff --git a/mex/sources/kalman/qt/f90/TU.f90 b/mex/sources/kalman/qt/f90/TU.f90
new file mode 100644
index 000000000..1bcdd253a
--- /dev/null
+++ b/mex/sources/kalman/qt/f90/TU.f90
@@ -0,0 +1,31 @@
+subroutine TU(X,T,U,n)
+!     COMPUTATIONAL SUBROUTINE: TU=T*U; T upper triangular matrix; U strictly upper triangular
+
+
+implicit none
+
+integer(4), intent(in) :: n
+real(8), intent(in) :: T(n,n)
+real(8), intent(in) :: U(n,n)
+
+real(8), intent(out) :: X(n,n)
+integer(4) :: i,j,k
+
+real(8) :: stemp
+
+
+X=0.0*T
+
+
+do i=1,n-1
+do j=i+1,n
+stemp = 0.0
+do  k = i,j-1
+stemp = stemp + T(i,k) * U(k,j)
+end do
+X(i,j)=stemp
+end do
+end do
+
+return
+end subroutine TU
diff --git a/mex/sources/kalman/qt/f90/TUt.f90 b/mex/sources/kalman/qt/f90/TUt.f90
new file mode 100644
index 000000000..050787044
--- /dev/null
+++ b/mex/sources/kalman/qt/f90/TUt.f90
@@ -0,0 +1,29 @@
+subroutine TUt(X,T,Ut,n)
+!     COMPUTATIONAL SUBROUTINE TUt=T*Ut T upper triangular; U striclty lower triangular
+
+implicit none
+
+integer(4), intent(in) :: n
+real(8), intent(in) :: T(n,n)
+real(8), intent(in) :: Ut(n,n)
+
+real(8), intent(out)  :: X(n,n)
+integer(4) :: i,j,k,h
+
+real(8) :: stemp
+
+
+X=0.0*T
+do i=1,n
+do j=1,n-1
+h=max(i,j)
+stemp = 0.0
+do  k = 0,n-h
+stemp = stemp + T(i,k+h) * Ut(k+h,j)
+end do
+X(i,j)=stemp
+end do
+end do
+
+return
+end subroutine TUt
\ No newline at end of file
diff --git a/mex/sources/kalman/qt/f90/TV.f90 b/mex/sources/kalman/qt/f90/TV.f90
new file mode 100644
index 000000000..bdc651f7e
--- /dev/null
+++ b/mex/sources/kalman/qt/f90/TV.f90
@@ -0,0 +1,27 @@
+subroutine TV(X,T,V,n)
+!     COMPUTATIONAL SUBROUTINE: TV=T*V T upper triangular; V vector
+
+implicit none
+
+integer(4), intent(in) :: n
+real(8), intent(in) :: T(n,n)
+real(8), intent(in) :: V(n,1)
+
+real(8), intent(out) :: X(n,1) 
+
+integer(4) :: i,k
+
+real(8) :: stemp
+
+
+
+do i=1,n
+stemp = 0.0
+do  k = 0,n-i
+stemp = stemp + T(i,k+i) * V(k+i,1)
+end do
+X(i,1)=stemp
+end do
+
+return
+end subroutine TV
diff --git a/mex/sources/kalman/qt/test/Makefile b/mex/sources/kalman/qt/test/Makefile
new file mode 100644
index 000000000..8e82c01db
--- /dev/null
+++ b/mex/sources/kalman/qt/test/Makefile
@@ -0,0 +1,65 @@
+# $Id: Makefile 531 2005-11-30 13:49:48Z kamenik $
+# Copyright 2005, Ondra Kamenik
+
+DEBUG = yes
+
+#LD_LIBS := -llapack -lcblas -lf77blas -latlas -lg2c
+
+CC_FLAGS :=  -DMATLAB -DWINDOWS  -DNO_BLAS_H -DNO_LAPACK_H \
+	-Wall -I../cc -I../sylv/cc -I../cc \
+	-Ic:/"Program Files"/MATLAB_SV71/extern/include  #-pg
+
+ifeq ($(DEBUG),yes)
+	CC_FLAGS := -DDEBUG $(CC_FLAGS) -g
+#	CC_FLAGS := -DTIMING_LOOP -DDEBUG $(CC_FLAGS) -g  
+	KALMANLIB := kalmanlib_dbg.a
+else
+#	CC_FLAGS := $(CC_FLAGS) -O2
+	CC_FLAGS := -DTIMING_LOOP $(CC_FLAGS) -O2
+	KALMANLIB := kalmanlib.a
+endif
+
+# Added by GP
+#    LDFLAGS := -llapack -lcblas -lf77blas -latlas -lg2c -lstdc++ -lmingw32
+        #LDFLAGS := -Wl,--library-path  $(LD_LIBRARY_PATH) 
+
+        LD_LIBS := -Wl,--library-path  \
+              -Wl,-L'f:/MinGW/lib' \
+        -Wl,-L"c:/Program Files"/MATLAB_SV71/extern/lib/win32/microsoft/ \
+        -Wl,-llibmex -Wl,-llibmx  -Wl,-llibmwlapack -Wl,-llibdflapack \
+        -lf95 -lg2c -lmingw32  -lstdc++ $(LDFLAGS) \
+        -Wl,-L'C:/MinGW/lib/gcc-lib/i686-pc-mingw32/4.0.4' -Wl,-L'C:/MinGW/lib' 
+
+#                  -Wl,-L'f:/CygWin/usr/local/atlas/lib' 
+#                  -Wl,-L'f:/CygWin/lib' 
+#                 $(LDFLAGS)
+#    LD_LIBS :=$(LDFLAGS)
+# end add
+
+#matrix_interface := GeneralMatrix Vector SylvException 
+#matobjs := $(patsubst %, ../sylv/cc/%.o, $(matrix_interface))
+#mathsource := $(patsubst %,  ../sylv/cc/%.h, $(matrix_interface))
+#matcppsource := $(patsubst %,  ../sylv/cc/%.cpp, $(matrix_interface))
+qtf90source := $(wildcard ../f90/*.f90)
+qtobjs := $(patsubst %.f90,%.o,$(qtf90source))
+cppsource := $(wildcard *.cpp)
+hsource := $(wildcard *.h)
+objects := $(patsubst %.cpp,%.o,$(cppsource))
+
+dummy.ch:
+	touch dummy.ch
+
+
+%.o: %.cpp $(hsource) $(cppsource)
+	c++ $(CC_FLAGS) -c $*.cpp
+
+qtamvm_exe.exe: qtamvm_exe.o $(qtobjs) $(hsource) $(cppsource) 
+	gcc $(CC_FLAGS) -o qtamvm_exe.exe qtamvm_exe.o ascii_array.o \
+	$(qtobjs) $(LD_LIBS) 
+
+all: $(objects) qtamvm_exe.exe # $(cppsource) $(hsource) $(kalmanhsource) $(kalmancppsource)
+
+clear:
+	rm -f *.o
+	rm -f *.a
+	rm -f *.{exe,dll}
diff --git a/mex/sources/kalman/qt/test/QT_tab b/mex/sources/kalman/qt/test/QT_tab
new file mode 100644
index 000000000..d5eeaa432
--- /dev/null
+++ b/mex/sources/kalman/qt/test/QT_tab
@@ -0,0 +1,11 @@
+0.153060262389052	0.001380442914015	0.001695719778633	-0.096495788866461	-0.006326378534057	0.000968937494692	0.00940481674898	-0.003274924112445	0.001481132666149	0.004112932319104	-0.012701307804641
+0	0.073952133671219	-0.001504207761519	0.057433281320695	0.009662045343528	0.006617685912409	0.013462091408161	0.007891283166287	0.00060603779672	0.003369481465008	0.018748481144843
+0	0	0.016722949277241	-0.043158612323585	-0.004579527748968	0.014353708194285	-0.002590188137874	0.015728887447029	-0.000497609190313	-0.000142760660804	-0.019083813953199
+0	0	0	0.265453840405543	0.00050417813553	0.002313090008356	0.027009998927894	-0.00302237841155	-0.006605759131604	-0.009638158849163	0.01228621141328
+0	0	0	0	0.040159616212026	-0.035581727103737	-0.016954216473738	-0.007949691886576	0.003628861926114	-0.001291923621455	0.039842615500298
+0	0	0	0	0	0.10084154296869	-0.001042311199417	-0.005382521361517	0.001696402254541	-0.000540820717391	0.023206920949216
+0	0	0	0	0	0	0.054131830294643	-0.020987700402038	0.000402689298563	-0.000018817482395	0.001917827090534
+0	0	0	0	0	0	0	0.14427448853052	-0.000904268210001	-0.004155316544672	0.018404810011954
+0	0	0	0	0	0	0	0	0.009309814517978	0.000447330038827	0.00529992488884
+0	0	0	0	0	0	0	0	0	0.005213746475416	0.012395564932675
+0	0	0	0	0	0	0	0	0	0	0.026158448255444
\ No newline at end of file
diff --git a/mex/sources/kalman/qt/test/QTmatlabTest.m b/mex/sources/kalman/qt/test/QTmatlabTest.m
new file mode 100644
index 000000000..3ae646523
--- /dev/null
+++ b/mex/sources/kalman/qt/test/QTmatlabTest.m
@@ -0,0 +1,45 @@
+size=100;
+AA=rand(size);
+BB=rand(size);
+[QT1 QT2 U1 U2]=qz(AA,BB);
+
+a=rand(size,1);
+
+Loops=10000
+
+% Calling QT without use of Sylv Vector and General Matrix
+t = clock;  
+[AAAA]=qtamvm(QT1,a,Loops);
+QTcpp_noSylv_TaInnerLoop_time=etime(clock, t)
+
+% Calling QT using of Sylv Vector and General Matrix
+t = clock;  
+[AAAA]=qtmvm(QT1,a,Loops);
+QTcppTaInnerLoop_time=etime(clock, t)
+
+t = clock;  
+[AAAA]=qtmvm_sub(QT1,a,Loops);
+QTcppTaInnerLoop_time=etime(clock, t)
+
+t = clock;  
+for tt=1:Loops%0
+[AAAA]=qtmvm_sub(QT1,a,1);
+end
+QTcppTaOuterLoop_time=etime(clock, t)
+
+t = clock;  
+[AAAA]=gmvm(QT1,a,Loops);
+GMcppTaInnrLoop_time=etime(clock, t)
+
+t = clock;  
+for tt=1:Loops%0
+[AAAA]=gmvm(QT1,a,1);
+end
+GMcppTaOuterLoop_time=etime(clock, t)
+
+t = clock;  
+for tt=1:Loops%0
+MTA=QT1*a;
+end
+matlabTa_time=etime(clock, t)
+
diff --git a/mex/sources/kalman/qt/test/aa_dat b/mex/sources/kalman/qt/test/aa_dat
new file mode 100644
index 000000000..cdf5c02e4
--- /dev/null
+++ b/mex/sources/kalman/qt/test/aa_dat
@@ -0,0 +1,11 @@
+0.130845306698951
+0.784832609671869
+0.213653301284573
+0.81536432324564
+0.021265975953734
+0.226338606586172
+0.447782600002914
+0.43948408700417
+0.687965758422965
+0.34628969510696
+0.996770191062048
\ No newline at end of file
diff --git a/mex/sources/kalman/qt/test/ascii_array.cpp b/mex/sources/kalman/qt/test/ascii_array.cpp
new file mode 100644
index 000000000..05c34a2d1
--- /dev/null
+++ b/mex/sources/kalman/qt/test/ascii_array.cpp
@@ -0,0 +1,84 @@
+// ascii_matrix.cpp 
+// Based on work of 2005, Ondra Kamenik
+
+#include "ascii_array.h"
+
+#include <stdio.h>
+#include <string.h>
+#include <iostream>
+
+#include <fstream>
+#include <string>
+
+// if the file doesn't exist, the number array is empty
+void 
+AsciiNumberArray::GetMX(const char* fname, int INrows, int INcols)
+  {
+  rows = 0;
+  cols = INcols;
+  
+  std::ifstream file(fname);
+  std::string line;
+  
+  data=(double*)calloc(INcols*INrows,sizeof(double));
+  while (getline(file, line)) 
+    {
+    rows++;
+    int icols = 0;
+    const char delims[] = " \r\n\t";
+    char* lineptr = strdup(line.c_str());
+    char* tok = strtok(lineptr, delims);
+    while (tok) 
+      {
+      icols++;
+      double item;
+      if (1 != sscanf(tok, "%lf", &item)) 
+        {
+        fprintf(stderr, "Couldn't parse a token %s as double.\n", tok);
+        exit(1);
+        }
+      data[(rows-1)*INcols+icols-1]=item;
+      tok = strtok(NULL, delims);
+      }
+    free(lineptr);
+    if (cols) 
+      {
+      if (cols != icols) 
+        {
+        fprintf(stderr, "Asserted a different number of columns.\n");
+        exit(1);
+        }
+      } 
+    else 
+      {
+      cols = icols;
+      }
+    }
+  }
+
+
+void 
+AsciiNumberArray::WriteMX()
+  {
+  std::ofstream outFile(strcat(fname, "_out"));
+  for (int i = 0; i < rows; i++)
+    {
+    for (int j = 0; j < cols; j++)
+      {
+      outFile << data[i*cols+j] << "\t";
+      }
+    outFile << std::endl;
+    }
+  outFile.close();
+  }
+
+
+void WriteMX(char* fname, double* data, int rows, int cols)
+  {
+  AsciiNumberArray OutArray;
+  OutArray.fname=fname;
+  OutArray.rows=rows;
+  OutArray.cols=cols;
+  OutArray.data=data;
+  OutArray.WriteMX();
+  }
diff --git a/mex/sources/kalman/qt/test/ascii_array.h b/mex/sources/kalman/qt/test/ascii_array.h
new file mode 100644
index 000000000..d91ad23bd
--- /dev/null
+++ b/mex/sources/kalman/qt/test/ascii_array.h
@@ -0,0 +1,15 @@
+// ascii_matrix.h 
+// Based on work of  Ondra Kamenik
+
+
+struct AsciiNumberArray {
+	int rows;
+	int cols;
+  double * data;
+	char* fname;
+  void GetMX(const char* fname,  int rows, int cols);
+  void WriteMX();
+};
+
+void WriteMX(char* fname, double* data, int rows, int cols);
+
diff --git a/mex/sources/kalman/qt/test/qtamvm_exe.cpp b/mex/sources/kalman/qt/test/qtamvm_exe.cpp
new file mode 100644
index 000000000..c407f9dbb
--- /dev/null
+++ b/mex/sources/kalman/qt/test/qtamvm_exe.cpp
@@ -0,0 +1,107 @@
+/*
+* Copyright (C) 2008-2009 Dynare Team
+*
+* This file is part of Dynare.
+*
+* Dynare is free software: you can redistribute it and/or modify
+* it under the terms of the GNU General Public License as published by
+* the Free Software Foundation, either version 3 of the License, or
+* (at your option) any later version.
+*
+* Dynare is distributed in the hope that it will be useful,
+* but WITHOUT ANY WARRANTY; without even the implied warranty of
+* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+* GNU General Public License for more details.
+*
+* You should have received a copy of the GNU General Public License
+* along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
+*/
+
+/****************************************************** 
+%
+% This provides an interface to QT f90 library by Andrea Pagano 
+% to multiply Quasi trinagular matrix (T) with a vector a
+%
+% function [a] = qtamvm(QT,a)
+%
+%  1. T1 = QT2T(QT;n) and Ld = QT2Ld(QT;n);
+%  2. Ta = LdV(Ld;a;n)+TV(T1;a;n).
+%
+% INPUTS
+%    T                      [double]    mm*mm transition matrix of the state equation.
+%    a                      [double]    mm state vector.
+%
+% OUTPUTS
+%    Tinverse               [double]    mm*mm transition matrix of the state equation.
+**********************************************************/
+
+#include "qt.h"
+#include <stdio.h>
+#include <math.h>
+#include "ascii_array.h"
+#include <iostream>
+#include <stdexcept>
+#include <malloc.h>
+
+int main(int argc, char* argv[])
+    {
+    
+    if (argc < 3 )
+      {
+      printf("Must have  min 2 input parameters.\n");
+      exit(1);
+      }
+  
+    try 
+      {
+      // make input matrices
+      int n=atoi(argv[3]);
+      double *T1, *Ld, *TV, * Ta ;
+      AsciiNumberArray QT, a;
+      QT.GetMX(argv[1],n,n);
+      a.GetMX(argv[2],n,1);
+      T1=(double *)calloc(n*n, sizeof(double));
+      Ld=(double *)calloc(n*n,sizeof(double));
+      TV=(double *)calloc(n, sizeof(double));
+      // create output and upload output data
+      Ta=(double *)calloc(n, sizeof(double));
+
+      
+#ifdef TIMING_LOOP
+      int loops=1;//000;
+      if (argc >3 )
+        loops = atoi(argv[2]);
+      for (int tt=0;tt<loops;++tt)
+        {
+#endif
+#ifdef DEBUG	
+      //  QT.print();
+#endif
+        // 1. T1 = QT2T(QT;n) and Ld = QT2Ld(QT;n);
+        qt2t_(T1, QT.data ,&n) ;
+        qt2ld_(Ld , QT.data,&n);
+        // 2. Ta = LdV(Ld;a;n)+TV(T1;a;n).
+        ldv_(Ta, Ld,a.data ,&n);
+        tv_(TV, T1 ,a.data,&n);
+//        Ta.add(1.0,TV);
+          for (int j=0; j<n;++j)
+            Ta[j]+=TV[j];
+
+
+#ifdef TIMING_LOOP
+        }
+      printf("QT array mvm: finished: %d loops\n",loops);
+#endif
+      // create output and upload output data
+      WriteMX(argv[2], Ta,n,1);
+      free(T1);
+      free(Ld);
+      free(TV);
+      free(Ta);
+      } 
+    catch (std::exception e) 
+      {
+      std::cout <<"Error" << std::endl;
+      }
+    
+  }; //main