diff --git a/filter.js b/filter.js new file mode 100644 index 0000000..7e81596 --- /dev/null +++ b/filter.js @@ -0,0 +1,135 @@ +var yb = [ 0.00010501, 0.00042004, 0.00063006, 0.00042004, 0.00010501]; +var ya = [ 1.00000000, -3.33787709, 4.20480923, -2.36821297, 0.50296098]; +var y_zi = [ 0.99989499, -2.33840213, 1.86577704, -0.50285597, 0.00000000]; + +var cb = [ 0.70911985, -2.83647938, 4.25471907, -2.83647938, 0.70911985]; +var ca = [ 1.00000000, -3.3302056 , 4.17977946, -2.34191808, 0.49401438]; +var c_zi = [-0.70911985, 2.12735954, -2.12735954, 0.70911985, 0.00000000]; + +var ib = [ 0.00605349, 0.01210697, 0.00605349]; +var ia = [ 1.00000000, -1.76816293, 0.79237688]; +var i_zi = [ 0.99394651, -0.78632339, 0.00000000]; + +var qb = [ 0.00278987, 0.00557974, 0.00278987]; +var qa = [ 1.00000000, -1.84512919, 0.85628867]; +var q_zi = [ 0.99721013, -0.85349880, 0.00000000]; + +function add_index_range(indices, beg, end) +{ + for (var i = beg; i <= end; ++i) + indices.push(i); + return indices; +} + +function add_index_const(indices, value, numel) +{ + while (numel--) + indices.push(value); + return indices; +} + +function subvector_reverse(vec, idx_end, idx_start) +{ + return vec.slice(idx_start, idx_end+1).reverse(); +} + +inline int max_val(const vectori& vec) +{ + return std::max_element(vec.begin(), vec.end())[0]; +} + +function filtfilt(B, A, X, zi) +{ + var len = X.length; // length of input + var na = A.length; + var nb = B.length; + var nfilt = (nb > na) ? nb : na; + var nfact = 3 * (nfilt - 1); // length of edge transients + + if (len <= nfact) { + console.log("Input data too short! Data must have length more than 3 times filter order."); + return X; + } + + // set up filter's initial conditions to remove DC offset problems at the + // beginning and end of the sequence + + var rows, cols; + //rows = [1:nfilt-1 2:nfilt-1 1:nfilt-2]; + rows = add_index_range([], 0, nfilt - 2); + if (nfilt > 2) + { + rows = add_index_range(rows, 1, nfilt - 2); + rows = add_index_range(rows, 0, nfilt - 3); + } + //cols = [ones(1,nfilt-1) 2:nfilt-1 2:nfilt-1]; + cols = add_index_const([], 0, nfilt - 1); + if (nfilt > 2) + { + cols = add_index_range(cols, 1, nfilt - 2); + cols = add_index_range(cols, 1, nfilt - 2); + } + // data = [1+a(2) a(3:nfilt) ones(1,nfilt-2) -ones(1,nfilt-2)]; + + var klen = rows.length; + var data = new Array(klen); + + data[0] = 1 + A[1]; var i, j = 1; + if (nfilt > 2) + { + for (i = 2; i < nfilt; i++) + data[j++] = A[i]; + for (i = 0; i < nfilt - 2; i++) + data[j++] = 1.0; + for (i = 0; i < nfilt - 2; i++) + data[j++] = -1.0; + } + + var leftpad = subvector_reverse(X, nfact, 1); + var _2x0 = 2 * X[0]; + for (i = 0; i < leftpad.length; ++i) {leftpad[i] = _2x0 - leftpad[i];} + + var rightpad = subvector_reverse(X, len - 2, len - nfact - 1); + var _2xl = 2 * X[len-1]; + for (i = 0; i < rightpad.length; ++i) {rightpad[i] = _2xl - rightpad[i];} + + var y0; + + var signal1 = leftpad.concat(X, rightpad); + + // Do the forward and backward filtering + y0 = signal1[0]; + zzi = zi.slice(0); // clone the array; + for (i = 0; i < zzi.length - 1; ++i) {zzi[i] = y0 * zzi[i];} + var signal2 = filter(B, A, signal1, zzi).reverse(); + zzi = zi.slice(0); + y0 = signal2[0]; + for (i = 0; i < zzi.length - 1; ++i) {zzi[i] = y0 * zzi[i];} + signal1 = filter(B, A, signal2, zzi); + + return subvector_reverse(signal1, signal1.length - nfact - 1, nfact); +} + +function filter(B, A, X, zi) +{ + + var input_size = X.length; + var filter_order = A.length + + // Initialize Y, the return variable + var Y = new Array(input_size); + for (var i = 0; i < input_size; ++i) {Y[i] = 0.0;} + + for (i = 0; i < input_size; ++i) { + var order = filter_order - 1; + while(order) { + if (i >= order) { + zi[order - 1] = B[order] * X[i - order] - A[order] * Y[i - order] + zi[order]; + } + --order; + } + Y[i] = B[0] * X[i] + zi[0]; + } + + return Y; +} diff --git a/lo-res-raw.png b/lo-res-raw.png new file mode 100644 index 0000000..1fe60e1 Binary files /dev/null and b/lo-res-raw.png differ diff --git a/test.html b/test.html new file mode 100644 index 0000000..c77cfdf --- /dev/null +++ b/test.html @@ -0,0 +1,199 @@ + + +
+r&&(.2>d&&(b[0].x+=1),.2>a&&(b[1].x+=1),.2>q&&(b[2].x+=1));n=0;for(p=this.vertices.length;n
c.y?this.quaternion.set(1,0,0,0):(a.set(c.z,0,-c.x).normalize(),b=Math.acos(c.y),this.quaternion.setFromAxisAngle(a,b))}}();
+THREE.ArrowHelper.prototype.setLength=function(a,b,c){void 0===b&&(b=.2*a);void 0===c&&(c=.2*b);this.line.scale.set(1,a,1);this.line.updateMatrix();this.cone.scale.set(c,b,c);this.cone.position.y=a;this.cone.updateMatrix()};THREE.ArrowHelper.prototype.setColor=function(a){this.line.material.color.set(a);this.cone.material.color.set(a)};
+THREE.BoxHelper=function(a){var b=new THREE.BufferGeometry;b.addAttribute("position",new THREE.BufferAttribute(new Float32Array(72),3));THREE.Line.call(this,b,new THREE.LineBasicMaterial({color:16776960}),THREE.LinePieces);void 0!==a&&this.update(a)};THREE.BoxHelper.prototype=Object.create(THREE.Line.prototype);
+THREE.BoxHelper.prototype.update=function(a){var b=a.geometry;null===b.boundingBox&&b.computeBoundingBox();var c=b.boundingBox.min,b=b.boundingBox.max,d=this.geometry.attributes.position.array;d[0]=b.x;d[1]=b.y;d[2]=b.z;d[3]=c.x;d[4]=b.y;d[5]=b.z;d[6]=c.x;d[7]=b.y;d[8]=b.z;d[9]=c.x;d[10]=c.y;d[11]=b.z;d[12]=c.x;d[13]=c.y;d[14]=b.z;d[15]=b.x;d[16]=c.y;d[17]=b.z;d[18]=b.x;d[19]=c.y;d[20]=b.z;d[21]=b.x;d[22]=b.y;d[23]=b.z;d[24]=b.x;d[25]=b.y;d[26]=c.z;d[27]=c.x;d[28]=b.y;d[29]=c.z;d[30]=c.x;d[31]=b.y;
+d[32]=c.z;d[33]=c.x;d[34]=c.y;d[35]=c.z;d[36]=c.x;d[37]=c.y;d[38]=c.z;d[39]=b.x;d[40]=c.y;d[41]=c.z;d[42]=b.x;d[43]=c.y;d[44]=c.z;d[45]=b.x;d[46]=b.y;d[47]=c.z;d[48]=b.x;d[49]=b.y;d[50]=b.z;d[51]=b.x;d[52]=b.y;d[53]=c.z;d[54]=c.x;d[55]=b.y;d[56]=b.z;d[57]=c.x;d[58]=b.y;d[59]=c.z;d[60]=c.x;d[61]=c.y;d[62]=b.z;d[63]=c.x;d[64]=c.y;d[65]=c.z;d[66]=b.x;d[67]=c.y;d[68]=b.z;d[69]=b.x;d[70]=c.y;d[71]=c.z;this.geometry.attributes.position.needsUpdate=!0;this.geometry.computeBoundingSphere();this.matrix=a.matrixWorld;
+this.matrixAutoUpdate=!1};THREE.BoundingBoxHelper=function(a,b){var c=void 0!==b?b:8947848;this.object=a;this.box=new THREE.Box3;THREE.Mesh.call(this,new THREE.BoxGeometry(1,1,1),new THREE.MeshBasicMaterial({color:c,wireframe:!0}))};THREE.BoundingBoxHelper.prototype=Object.create(THREE.Mesh.prototype);THREE.BoundingBoxHelper.prototype.update=function(){this.box.setFromObject(this.object);this.box.size(this.scale);this.box.center(this.position)};
+THREE.CameraHelper=function(a){function b(a,b,d){c(a,d);c(b,d)}function c(a,b){d.vertices.push(new THREE.Vector3);d.colors.push(new THREE.Color(b));void 0===f[a]&&(f[a]=[]);f[a].push(d.vertices.length-1)}var d=new THREE.Geometry,e=new THREE.LineBasicMaterial({color:16777215,vertexColors:THREE.FaceColors}),f={};b("n1","n2",16755200);b("n2","n4",16755200);b("n4","n3",16755200);b("n3","n1",16755200);b("f1","f2",16755200);b("f2","f4",16755200);b("f4","f3",16755200);b("f3","f1",16755200);b("n1","f1",16755200);
+b("n2","f2",16755200);b("n3","f3",16755200);b("n4","f4",16755200);b("p","n1",16711680);b("p","n2",16711680);b("p","n3",16711680);b("p","n4",16711680);b("u1","u2",43775);b("u2","u3",43775);b("u3","u1",43775);b("c","t",16777215);b("p","c",3355443);b("cn1","cn2",3355443);b("cn3","cn4",3355443);b("cf1","cf2",3355443);b("cf3","cf4",3355443);THREE.Line.call(this,d,e,THREE.LinePieces);this.camera=a;this.matrix=a.matrixWorld;this.matrixAutoUpdate=!1;this.pointMap=f;this.update()};
+THREE.CameraHelper.prototype=Object.create(THREE.Line.prototype);
+THREE.CameraHelper.prototype.update=function(){var a,b,c=new THREE.Vector3,d=new THREE.Camera,e=function(e,g,h,k){c.set(g,h,k).unproject(d);e=b[e];if(void 0!==e)for(g=0,h=e.length;g