Some coming updates

This will be just a short post to show a couple of the new graph types that will be available in the next release of OT/JSvg. Both of them are the result of a project I've been working on for the orthopedics guys at duPont. The first type, a bar graph with error bars, is rather standard, but was not possible using OT/JSvg without all your own extra coding.

It's nothing fancy, but the syntax is rather straight forward:

      drawBarSeriesErrorPlot(xnames, yarr,yerr,fn);

In this case, xnames is a 1d array of labels for the x-axis, yarr is a 2d array of y-values for plotting in series form, yerr is the corresponding 2d array of error values to be added or subtracted from yarr, and fn is the filename to save.

But moving along, we have something, which is to my mind quite a lot more impressive: mapping! Using shapefiles converted from the ogr2ogr utility (used by GRASS and other GIS systems), demographics files from the US census, and then event files (for example number of crimes or accidents), one can then paint a graph as follows.

This map doesn't include events at the current time as that data is not yet public. If you know Delaware, you'll also note that parts of New Castle are missing. That's because the MULTIPOLY parsing isn't working exactly right yet. Other features are an ability to overlay other layers for visualization. So you can draw the map, color the map dependent on the feature, insert points for events, and insert overlay data (gradients, heat maps, etc) to show additional trends.

     drawMap(xCoords, yCoords, refLongs, refLats, overlayData, overlayCT, xPoints, yPoints, coordType, fn, legendStuff) {
    // xCoords - x values or longitude data to be converted (2D)
    // yCoords - y values or latitude data to be converted (2D)
    //      For each series above, a line color and fill type are needed
    // refLong - a reference longitude for drawing map.  Greenwich is far away and 
    //      causes distortion of map over great distances
    // overlayData - 3D array containing data for coloration
    //      This will always be a rectangular array whose upper left position
    //      corresponds to the minimum xCoord and minimum yCoord, and whose
    //      lower right position corresponds to the maximum xCoord and yCoord.
    //      This routine will do the registration of overlayData
    // //overlayCorners - 3D array containing coord data for overlayData
    //      =[ [ [xTop1,xBot1],[yTop1,yBot1] ], ...];
    //      These should be within the xCoords/yCoords range
    // overlayCT - colortable for overlayData 
    // xPoints - x values of points to overlay on top of everything else (2D)
    // yPoints - y values of points to overlay on top of everything else (2D)
    // coordType - 'grid' values or 'longlat' data  if the latter, then convert
    //      For xCoords, yCoords, and xPoints, yPoints, all should be of same type
    // fn - filename for svg file.

For you GIS types, this might not be that impressive, but what comes along with this is a regional ANOVA analysis, Getis-Ord* statistics by point and region, and some frequency analyses by region. It's not quite ready for prime time, at this point, but as I settle on the API, this function and its relatives will make it into the main trunk.

More shortly as there is time.

Cheers!

Snowday data modeling

It's a snowday here in South-Eastern Pennsylvania, so I'm off from work where I teach physics, astronomy and Taijiquan. As the weather is so unpredictable in these parts in the spring semester, we've come up with some computationally based labs to help students better understand the various physics phenomena they've been studying in the event that labs are canceled. Before I go any further, let me say that we attempt to teach the students the basics of JavaScript and give them a number of exercises before going any further.

Today's exercise is a bit more complicated. We're going to generate electric field lines for a distribution of charges. There are a number of ways to do this, but the way we will pick is based on the kinematics of a test charge. In principle, test charges move directly along field lines, and so if one maps out their paths, one also gets a sense of the field line direction. Sounds pretty easy, right? Basically yes, but the devil is in the details. Here is the general algorithm:

1. Make arrays to contain the positions of the charge distribution. These charges are point charges that are fixed in place so we will only need to specify a list of x positions and y positions. Therefore, 1D arrays will due. Call these variables: xpos and ypos, respectively.
2. Make arrays to contain the positions of the test charges. These charges are point charges that move. If we want to draw them, we need to have their positions at a number of different times. Therefore, these will be 2D arrays. The first dimension will refer to the test charge, and the second dimension will refer to the time. Call these variables: x and y.
3. Pick a reasonable time increment dti, and test charge mass masst.
4. Start with initial velocities v0x = 0; and v0y = 0;.
5. Calculate the x and y components of force on a given test charge due to ALL the charges in the charge distribution. Call these values Fx and Fy. These are determined by the Coulomb force.
6. Calculate the accelerations of each of the given test charges based on their masses (all masst. The accelerations will be ax = Fx/masst and ay = Fy/masst.
7. Update the positions of each of the test charges using the regular equations of kinematics. The x and y positions are given by
8. x[charge][time] = x[charge][time - 1] + v0x*dti + 0.5*ax*dti*dti
9. y[charge][time] = y[charge][time - 1] + v0y*dti + 0.5*ay*dti*dti
Remember: we only need to store the x and y positions to draw the field lines.
10. Update the x and y velocities of each test charge:
11. v0x = v0x + ax*dti
12. v0y = v0y + ay*dti
13. Go back to step 5. Rinse, wash, repeat.

In code, we will start by layout out the initial conditions.

var q = 1.0e-6;     // C (magnitude of each charge in charge dist)
var qt = 1.0e-10;   // C (test charge magnitude)
var masst = 1.0e-10;// kg (test charge mass)
var k = 8.99e9;     // N*m^2/C^2  (Coulomb's constant)

// Positions of charge array
var ypos = new Array();
var xpos = new Array();

var N = 3;         // number of charges in charge distrib
var scale = 100;   // used in setting test charge positioning

var extra = 4      // how many MORE test charges than charges in distrib
var Nt = 20;       // Number of time increments
var dt = 5.0e-6;   // seconds (time jump per increment) 

var Ntest = N+2*extra;  // how many test charges 
                        //   2*extra on top AND bottom

Next we set the charge distribution positions. In this case, it's going to be a vertical line at the origin, not quite centered in y.

// These set the charge positions
for (var i=0; i< N; i++) {
    ypos[i] = (i-N/2)/scale;
    xpos[i] = 0;           
}

Having set these up, we then create an array full of our test charges. We'll space them vertically along side our charge distribution, and let them be offset just a bit in x. Why the offset? Because the Coulomb force is infinite at a distance of zero. But before we go on, let's set them according to the scale of the fixed charge distribution. dy is the difference in height between any two charges in this distribution.

var dy = Math.abs(ypos[1]-ypos[2]);   // Difference in charge positions

Now, create the arrays for the test charges, and set an initial position for the lowest test charge. We'll set initial positions for the rest of them shortly.

// These arrays will be 2D and will store the positions of the test charges
// which are a proxy for the field lines.
var x = new Array();    
var y = new Array();

var x0 = 0.005;  // m
var y0 = ypos[0]-extra*dy; // m   

We finally come to the looping which will iterate over test charges and over the number of time increments we chose. First a note on Coulomb's Law. While it is typically written as a 1/r² law, if one applies vector math properly, it is possible to get components without using trigonometry explicitly as indicated in this sketch.

In that case, there is an r³ term in the denominator. The distance between a charge from the distribution and a test charge is given by

            // Distance (x and y) between field point and charge in array
            var deltaX = x[i][j-1] - xpos[p];
            var deltaY = y[i][j-1] - ypos[p];

The distance r³ is then given by this expression

            var r3 = Math.pow(deltaX*deltaX + deltaY*deltaY,1.5);

Then, following the rules of vector math, we (finally) get a form for the Coulomb force as a function of distance:

            //  Fx = k*q1*q2*x/r3 and
            //  Fy = k*q1*q2*y/r3
            Fx += k*qt*q*deltaX/r3;
            Fy += k*qt*q*deltaY/r3;

At this point, you're ready to see the full code for the positioning of the test charges at any time during the simulation. Note that there is also some code to draw arrows.

for (var i=0; i< Ntest; i++) {
    // Loop through Nt test charges, each starting from a point x0 to right of
    // the charge array but at the same height as the charge to the left.
    
    x[i] = new Array(); // Create 2nd dimension of the array
    y[i] = new Array();
    
    y0 += dy;   // m
    
    // Initial positions are known
    x[i][0] = x0;
    y[i][0] = y0;
    
    dti = dt;   //(Math.abs(Ntest/2-i)+0.01)*dt;    
    
    // The test charges start from rest
    var v0x = 0;    // m/s
    var v0y = 0;    // m/s

    // Start loop at j=1 since initial positions are known
    // This is a loop on time increments
    for (var j=1; j< Nt; j++) {
        
        // Now loop through all the charges to determine the Force in x and y
        var Fx = 0;
        var Fy = 0;
            
        // This is a loop on position for a given test charge
        for (var p=0; p< N; p++) {
            // Distance (x and y) between field point and charge in array
            var deltaX = x[i][j-1] - xpos[p];
            var deltaY = y[i][j-1] - ypos[p];
            // r3 = r cubed from Pythagorean theorem ( 3/2 power = 1.5)
            var r3 = Math.pow(deltaX*deltaX + deltaY*deltaY,1.5);
            // These have the form 
            //  Fx = k*q1*q2*x/r3 and
            //  Fy = k*q1*q2*y/r3
            Fx += k*qt*q*deltaX/r3;
            Fy += k*qt*q*deltaY/r3;
        }
        
        var ax = Fx/masst; // x-acceleration of test charge -changes
        var ay = Fy/masst; // y-acceleration of test charge -also changes
                
        x[i][j] = x[i][j-1] + v0x*dti + 0.5*ax*dti*dti;
        y[i][j] = y[i][j-1] + v0y*dti + 0.5*ay*dti*dti;
        
        v0x = v0x + ax*dti;
        v0y = v0y + ay*dti;
 
    }

    if (showArrows) {
        // drawArrowAtEnd(xp,yp,vx,vy,size)
        // x[i][j-1], y[i][j-1], ax, ay are last calculated values in preceding loop.    
        var tmp = drawArrowAtEnd(x[i][j-1],y[i][j-1],v0x,v0y,dy/20);
        var ind = 0;
        for (var jj = j; jj< j+3; jj++) {
            x[i][jj] = tmp[0][ind];
            y[i][jj] = tmp[1][ind];
            ind++;
        }
    }
    
    optsSeries[i].itype=1;        // Be sure it's a line rather than a point.
    optsSeries[i].thick=1;        // Make this 1px thick
    optsSeries[i].stroke="red";   // Make this path red
}

The last bit is code to figure out how to draw arrows at the end of the test charge paths. The inputs are the last calculated point for the test charge, the last velocities, and a length (dy/20) in this case. The output of drawArrowAtEnd() is a 2D array [xpts,ypts] which contains three x and three y points: the x and y coordinates of each side of the arrow, and the point. There is also code here for coloring the paths of the test charges. [Note that drawArrowAtEnd() is not part of OT/JSvg, but is included at the end of this article.

Having done all the calculations, we're now ready to plot. We're going to use our well known drawSeriesPlot() command, so we need to use a 2D array for both x and y. Our test charges are already configured like this, so let's see what happens if we just call our plot command.

Ugh! We certainly don't want a legend, and we can't see the charge distribution. Also, there are no labels. Let's do the easy stuff first and deal with these issues.

optsGraphTitles[0].name = 'x (m)';
optsGraphTitles[1].name = 'y (m)';
optsGraphTitles[2].name = 'Field lines with N = '+N+' charges';

window.optsLegend = false;

OK, this is better, but we also want to add the charge distribution. To do this, we'll just add them to the x and y arrays. Yes, these were originally for the test charges, but now that the calculations are done, we can just add in one more series.

// Also plot the charges

xl = x.length;
x[xl] = xpos;
y[xl] = ypos;

optsSeries[xl].itype=0;            // Make these points
optsSeries[xl].thick=2;            // Make them 2 wide
optsSeries[xl].stroke="blue";      // Color them blue (outline)
optsSeries[xl].color="blue";       // Make the fill blue also

Excellent! Now we are ready to play. Notice how the field diverges? What happens if you add more charges. Let's go with 20 charges, instead of 3. Set N = 3 at the top of the script.

Notice that the field lines in the middle have become more parallel. Want to make the arrows bigger? Change the final parameter of drawArrowAtEnd(). How about dy/2?

Here is the final script used to produce the final graph in its entirety:

/************ Only change things in this section ***********************/

var showArrows = true;

var q = 1.0e-6; // C
var qt = 1.0e-10;   // C (test charge magnitude)
var masst = 1.0e-10;    // kg (test charge mass)
var k = 8.99e9;     // N*m^2/C^2

// Position of charge array
var ypos = new Array();
var xpos = new Array();

var N = 20; // number of charges
var scale = 100;

var extra = 4
var Nt = 20;    // Number of time increments
var dt = 5.0e-6;    // seconds

var Ntest = N+2*extra;

/***********************************************************************/


function drawArrowAtEnd(xp,yp,vx,vy,size) {
    // xp and yp are the end point coordinates
    // vx and vy are vector components - used to get the angle
    // size is the size in x-y space of the arrow end points
    
    // The arrow is drawn by drawing points starting at the end point
    // and doing the lower "fletching", moving back to the end point and
    // then doing the upper "fletching".
    
    var r2d = 180/Math.PI;
    var d2r = 1/r2d;
    
    // calculate the slope of the vector
    if (xp != 0) {
        if (yp/xp < 1000) {
            // To ensure sanity
            var theta0 = r2d*Math.atan(yp/xp);
            var theta1 = (theta0-135)*d2r;
            var theta2 = (theta0+135)*d2r;
        } else { // do the same thing as if xp = 0
            if (yp < 0) {   // It points down
                theta1 = 135*d2r;
                theta2 = 45*d2r;
            } else {    // It points up
                theta1 = 45*d2r;
                theta2 = 135*d2r;
            }
        }

    } else {
        if (yp < 0) {   // It points down
            theta1 = 135*d2r;
            theta2 = 45*d2r;
        } else {    // It points up
            theta1 = 45*d2r;
            theta2 = 135*d2r;
        }    
    }
    console.log(theta1,theta2)
    var dx1 = size*Math.cos(theta1);
    var dy1 = size*Math.sin(theta1);
    var dx2 = size*Math.cos(theta2);
    var dy2 = size*Math.sin(theta2);

    var xout = [xp+dx1,xp,xp+dx2];
    var yout = [yp+dy1,yp,yp+dy2];
    return [xout,yout]; // 2d array containing arrow
}


// These set the charge positions
for (var i=0; i< N; i++) {
    ypos[i] = (i-N/2)/scale;
    xpos[i] = 0;           
}

var dy = Math.abs(ypos[1]-ypos[2]);   // Difference in charge positions

// These arrays will be 2D and will store the positions of the test charges
// which are a proxy for the field lines.
var x = new Array();    
var y = new Array();

var x0 = 0.005;  // m
var y0 = ypos[0]-extra*dy; // m   


// Now start a test charge at a distance of x0 from the array of charges and 
// see where it goes.  This will indicate the direction of the field line.  To
// do this, one must first calculate the net field at the point in order to get
// the force on the charge and thereby its acceleration.
// Use equation x(i+1) = x(i) + v0x*ttot + 0.5*ay(i)*ttot^2 (pseudo code)
// Ditto with y(i+1)
// In this case the particle is likely to go VERY fast.  Choose a very small dt:

for (var i=0; i< Ntest; i++) {
    // Loop through Nt test charges, each starting from a point x0 to right of
    // the charge array but at the same height as the charge to the left.
    
    x[i] = new Array(); // Create 2nd dimension of the array
    y[i] = new Array();
    
    y0 += dy;   // m
    
    // Initial positions are known
    x[i][0] = x0;
    y[i][0] = y0;
    
    dti = dt;   //(Math.abs(Ntest/2-i)+0.01)*dt;    
    
    // The test charges start from rest
    var v0x = 0;    // m/s
    var v0y = 0;    // m/s

    // Start loop at j=1 since initial positions are known
    // This is a loop on time increments
    for (var j=1; j< Nt; j++) {
        
        // Now loop through all the charges to determine the Force in x and y
        var Fx = 0;
        var Fy = 0;
            
        // This is a loop on position for a given test charge
        for (var p=0; p< N; p++) {
            // Distance (x and y) between field point and charge in array
            var deltaX = x[i][j-1] - xpos[p];
            var deltaY = y[i][j-1] - ypos[p];
            // r3 = r cubed from Pythagorean theorem ( 3/2 power = 1.5)
            var r3 = Math.pow(deltaX*deltaX + deltaY*deltaY,1.5);
            // These have the form 
            //  Fx = k*q1*q2*x/r3 and
            //  Fy = k*q1*q2*y/r3
            Fx += k*qt*q*deltaX/r3;
            Fy += k*qt*q*deltaY/r3;
        }
        
        var ax = Fx/masst; // x-acceleration of test charge -changes
        var ay = Fy/masst; // y-acceleration of test charge -also changes
                
        x[i][j] = x[i][j-1] + v0x*dti + 0.5*ax*dti*dti;
        y[i][j] = y[i][j-1] + v0y*dti + 0.5*ay*dti*dti;
        
        v0x = v0x + ax*dti;
        v0y = v0y + ay*dti;
 
    }

    if (showArrows) {
        // drawArrowAtEnd(xp,yp,vx,vy,size)
        // x[i][j-1], y[i][j-1], ax, ay are last calculated values in preceding loop.    
        var tmp = drawArrowAtEnd(x[i][j-1],y[i][j-1],v0x,v0y,dy/2);
        var ind = 0;
        for (var jj = j; jj< j+3; jj++) {
            x[i][jj] = tmp[0][ind];
            y[i][jj] = tmp[1][ind];
            ind++;
        }
    }
    
    optsSeries[i].itype=1;
    optsSeries[i].thick=1;
    optsSeries[i].stroke="red";
}

// Also plot the charges

xl = x.length;
x[xl] = xpos;
y[xl] = ypos;

optsSeries[xl].itype=0;
optsSeries[xl].thick=2;
optsSeries[xl].stroke="blue";
optsSeries[xl].color="blue";

optsGraphTitles[0].name = 'x (m)';
optsGraphTitles[1].name = 'y (m)';
optsGraphTitles[2].name = 'Field lines with N = '+N+' charges';

window.optsLegend = false;
/*
window.xprec = 4;
*/
// OK, let's plot this stuff

drawSeriesPlot(x,y,0,0,'fieldLines.svg');

Until next time: cheers!

More Fitting Examples

In the last post, I went through (in some detail) the details around how to do fits to data. In doing so, I included a fairly detailed discussion of the JavaScript array object. In this post, I'm going to be a bit less detailed and give a few examples of the use of the doFit() function. This is not, by any means, the only way to do fitting in JSvg, but it does provide a simple wrapper that hides the more complicated fitting methods from the user. Once the Analysis Console is up, if one types:

doFit.toSource()

the following code will appear. In fact, much more will appear, but I've shown the headers of the code here.

"function doFit(xArr,yArr,fitN,prec) {
  // This should return a fit array and the equations of fit, and coefficients
  // pass 1d xArr and 1d yArr
  // fitN = 1 - linear
  // 	2 - quadratic
  //	3 - cubic
  //	4 - sinusoidal"

In much of the JSvg codebase, the code is commented so as to provide the user with additional help as to how to use the various functions. In its current state, doFit() supports linear, quadratic, cubic, and sinusoidal fitting. (Note that in the next release of OT/JSvg, it will also support logarithmic, exponential, and power law fitting, but more on that in another post). So to use the function one needs a vector of x-data, a vector of y-data, a number specifying which fit to do, and a precision. If precision isn't specified, all output numbers will be shown to the third digit past the decimal place.

So, without further ado, here are a couple of examples. Presume that I have x-data and y-data as follows:

    var xdata = [1.3333,1.3381,1.3429,1.3476,1.3524,1.3571,1.3619,1.3667,1.3714,1.3762,1.3810,1.3857,1.3905,1.3952,1.4000,1.4048,
1.4095,1.4143,1.4190,1.4238,1.4286,1.4333,1.4381,1.4429,1.4476,1.4524,1.4571,1.4619,1.4667,1.4714,1.4762,1.4810,
1.4857,1.4905,1.4952,1.5000,1.5048,1.5095,1.5143,1.5190,1.5238,1.5286,1.5333,1.5381,1.5429,1.5476,1.5524,1.5571,
1.5619,1.5667,1.5714,1.5762,1.5810,1.5857,1.5905,1.5952,1.6000,1.6048,1.6095,1.6143,1.6190,1.6238,1.6286,1.6333,
1.6381,1.6429,1.6476,1.6524,1.6571,1.6619,1.6667,1.6714,1.6762,1.6810,1.6857,1.6905,1.6952,1.7000,1.7048,1.7095,
1.7143,1.7190,1.7238,1.7286,1.7333,1.7381,1.7429,1.7476,1.7524,1.7571,1.7619,1.7667,1.7714,1.7762,1.7810,1.7857,
1.7905,1.7952,1.8000,1.8048,1.8095,1.8143,1.8190,1.8238,1.8286,1.8333,1.8381,1.8429,1.8476,1.8524,1.8571,1.8619,
1.8667,1.8714,1.8762,1.8810,1.8857,1.8905,1.8952,1.9000,1.9048,1.9095,1.9143,1.9190,1.9238,1.9286,1.9333,1.9381,
1.9429,1.9476,1.9524,1.9571,1.9619,1.9667,1.9714,1.9762,1.9810];

    var ydata = [0.986,1.003,1.015,1.027,1.04,1.052,1.06,1.073,1.085,1.093,1.105,1.114,1.126,1.134,1.142,1.155,1.163,1.167,1.175,
1.184,1.192,1.2,1.208,1.216,1.225,1.229,1.237,1.245,1.249,1.253,1.262,1.266,1.27,1.278,1.282,1.286,1.29,1.295,1.299,
1.303,1.307,1.307,1.311,1.315,1.315,1.315,1.319,1.323,1.323,1.323,1.323,1.327,1.327,1.327,1.327,1.327,1.327,1.327,
1.327,1.327,1.323,1.323,1.323,1.319,1.315,1.315,1.315,1.311,1.307,1.307,1.303,1.299,1.295,1.29,1.286,1.282,1.278,
1.27,1.266,1.262,1.253,1.249,1.245,1.237,1.229,1.225,1.216,1.212,1.204,1.196,1.188,1.179,1.171,1.163,1.155,1.147,
1.138,1.13,1.118,1.11,1.097,1.089,1.077,1.064,1.044,1.044,1.036,1.023,1.011,0.999,0.986,0.974,0.962,0.949,0.933,
0.921,0.908,0.892,0.879,0.867,0.851,0.838,0.822,0.81,0.793,0.777,0.76,0.744,0.727,0.711,0.695,0.678,0.662,0.645,
0.625,0.608,0.592];

To plot the data with a fit, we must first construct a fit. We do it as follows (blindly assuming that these data follow a linear trend:

    var fitOrder = 1;
    var fitPrec = 3;
    fitObj1 = doFit(xdata,ydata,fitOrder,fitPrec);

Just because, at this point, we don't know that the data *are* actually linear, let's try fitting a quadratic as well:

    var fitOrder = 2;
    var fitPrec = 3;
    fitObj2 = doFit(xdata,ydata,fitOrder,fitPrec);

So now there are two "fitObjs" which contain fit data. What do they contain? Recall that doFit() returns the following values.

     
     //  fitObj[0] = 1D array containing the predicted y values for each x data point
     //  fitObj[1] = equation of fit (t is default variable)  (String)
     //  fitObj[2] = 1D array containing the coefficients of fit starting with highest order term 

So let's put the rest of our script together and plot.

    // Put x and y data and fits into output plotting variables:
    var xout = [xdata, xdata, xdata]; // xdata are the coordinate pairs for the fit data from doFit();
    var yout = [ydata, fitObj1[0], fitObj2[0]];  // matching y data

    // Now for graphics parameters:
    optsSeries[0].name = 'data';  // Call it what you will, but this is our original data
    optsSeries[1].name = fitObj1[1];  // Label for the first fit
    optsSeries[2].name = fitObj2[1];  // Label for the second fit

    // Fits look better with a line rather than points
    optsSeries[1].itype = 1;
    optsSeries[2].itype = 1;

    // And sometimes it's nice to adjust the colors
    optsSeries[1].stroke = 'blue';
    optsSeries[2].stroke = 'black';

    // How about some labels?
    optsGraphTitles[0].name = 'Time (s)';
    optsGraphTitles[1].name = 'y position (m)';
    optsGraphTitles[2].name = 'y vs t';
    optsGraphTitles[3].name = 'for a tossed ball';

    // Name the file, don't forget .svg Windows users!
    var fn = 'fallingMotion.svg';

    // Now plot it!
    drawSeriesPlot(xout,yout, 0,0, fn);

The resultant graph looks like this (after clicking the legend and title to adjust the positioning).

For the physics crowd, you will clearly see that these data originated from a ball drop and, as such, the quadratic fits much better. In fact these data result from the "Falling Motion" lab included in OpenTrack's tracking portion of the program.

A new release of OT/JSvg is due out shortly, and when it comes I will go over some of the new features.

Cheers!

Fitting a Line to Data with JSvg (and a brief discussion of JavaScript Arrays)

In the last post, we covered the basics of plotting a scattergraph using JSvg. In this post we will fit a line to the data from the last post as well as review some concepts.
When plotting anything on most 2D graphs (not including pie charts), one needs x-y coordinate pairs. This is even true in bar graphs and histograms (although generally the programs used to plot these other graph types will automatically generate the x-axis data). In JavaScript, the easiest way to store such data is in arrays. A one-dimensional JavaScript array can be likened to a list. It can be a list of numbers:

     var myArray = [1,2,3,4,5];

or

     var myArray = [1,5,23,43,2,7];

or even

     var myArray = [,,,2,3,4];

Where in the last case, the first three elements of the array are undefined. Arrays can also be lists of other things such as strings (representing characters or words):

     var myArray = ['Bob','Chris','a','b','c'];

In fact, arrays can be lists of other arrays:

     var arr1 = [1,2,3];
     var arr2 = [4,5,6];
     var arr3 = ['a,'b','c','d'];
     var myArray = [arr1,arr2,arr3];

Notice in the last case that the outer array (myArray) contains the other previously defined arrays. Also notice that not all arrays contain the same sort of thing. The first two inner arrays contain numbers, while the last one contains strings (in this case letters). Notice also, that the inner arrays are not all of the same length. This is allowable in JavaScript, and once you get used to it, it brings a great deal of flexibility to the table. In this last example, the array, myArray, could also be written like this:

     var myArray = [ [1,2,3], [4,5,6], ['a','b','c','d'] ];

or for clarity:

     var myArray = [ 
                     [1,2,3], 
                     [4,5,6], 
                     ['a','b','c','d'] 
                   ];

although this latter form takes up more space. To see how to access an element of this array, we could type into the Web Console:

     >> myArray = [ [1,2,3], [4,5,6], ['a','b','c','d'] ];

which would return:

     >> [objectArray]

Clicking on [objectArray] in the console results in the following:

At the right hand side, the contents of the [objectArray] are shown.  You can see that it consists of three inner (or sub) arrays and that its length is three.  The indices of the sub arrays are shown at the left  as 0, 1, and 2.  If we click on the third item in the array (index of 2), we see the following:

which exposes the content of the third inner array.  Let's say that we wish to access the content of the first inner array (index 0) from command line, we could type:

     myArray[0]

Which would result in another [object Array]. If we wished to access the first element of this array, we would type:

     myArray[0][0]

And in this case, since it is just a number, we would see its value (1) returned above in the console.

In general, in the future, we will use the language "zeroeth element" to refer to index 0 of an array, "first element" to refer to index 1 of an array and so on.  This is in keeping with JavaScript's convention to start indexing arrays at zero rather than one. 

After this rather lengthy preamble, we are ready to discuss fitting a line to the data from before.  The previous script was as follows:

     
    var R = 0.02275;   // m
    var Ncoils = 100;  // Number of coils
    var Ntrials = 6;   // Number of trials 
                   //   (presume we sample resistance every 100 coils)
    var pi = Math.PI;  // Use pi as a short hand for the JavaScript constant Math.PI
    var L = new Array();   // Define the array.  It's empty now.  
    for (var i=0; i< Ntrials; i++) {
         var Ntot = (i+1)*Ncoils;  // Total number of coils this sample.
         L[i] = 2*pi*R*Ntot;   // This is the ith value of length in the array
    }

    resist = [0.6,1.1,1.6,2.1,2.7,3.5];   // Ohms

    // Graphics commands
    optsGraphTitles[0].name = 'Length (m)';    // x-axis title
    optsGraphTitles[1].name = 'Resistance (\\Omega )';  // y-axis title
    optsGraphTitles[2].name = 'Resistance vs Length';   // main title
    optsGraphTitles[3].name = 'for a coil of wire';     // sub title

    optsLegend = false;

    drawSeriesPlot([L],[resist], 0,0, 'resistivityPlot.svg');

The idea is that we would be fitting a line to the x-y data represented by the array L and resist respectively. The big picture for physics students is that we are trying to determine the resistivity of the wire based on measurements of electrical resistance at given lengths, and knowledge of the wire's cross-sectional area (in this case: A = 4.64e-7 m²). The governing equation here is:

R = ρ L/A

where ρ is the resistivity. So if y corresponds to R and x corresponds to L, then ρ/A is the slope.
For quick fits, we can use the JSvg function "doFit()" which returns, you guessed it, an array! The syntax of the command is:

     
     fitObj = doFit(xvector, yvector, fitOrder, precisionOfoutput);

where xvector is the x-array, yvector is the y-array, fitOrder is 1 for linear, 2 for quadratic, 3 for cubic (and so on), and precisionOfoutput is the number of decimal places to show in the output. The function doFit returns an array that contains the following:

     
     fitObj[0] = 1D array containing the predicted y values for each x data point
     fitObj[1] = equation of fit (t is default variable)  (String)
     fitObj[2] = 1D array containing the coefficients of fit starting with highest order term 

There is a lower level function called by doFit() named polyCor() which contains much of the statistical information about the fit, but that's for another post.
So to modify the code from the prior and get a fit, we need to call doFit() with the appropriate arguments, and then extract the information we wish to display. Finally, we need to put it on our graph.

     
    var fitObj = doFit(L,resist,1,4);

We can then extract the objects we need and calculate the resistivity as follows:

     
// ********************* New code block ***********************************

    var fitObj = doFit(L,resist,1,4); 
    var fitData = fitObj[0];
    var fitEqu = fitObj[1];
    var slope = fitObj[2][0];

    var A = 4.64e-7;    // m^2
    var rho = slope*A;

    fitEqu = 'R = '+fitEqu.replace('t','L');
// ******************** End new code block ********************************

In the last equation, we have done some other formatting. By default the equation uses "t" as the x-variable. Let's replace it with "L", and then prepend "R = " on the front end to actually make it into an equation.
The only other item of substance is getting this information onto the graph. Recall that the function drawSeriesPlot() takes 2D arrays as its first two arguments. We recast this line as follows:

     
    drawSeriesPlot([L,L],[resist,fitData], 0,0, 'resistivityPlot2.svg');

What we have done is added another 1D array to both the x and y arrays in order to graph a second series on the plot. Below is the final code with all the bells and whistles.

     
    var R = 0.02275;   // m
    var Ncoils = 100;  // Number of coils
    var Ntrials = 6;   // Number of trials 

                       //   (presume we sample resistance every 100 coils)
    var pi = Math.PI;  // Use pi as a short hand for the JavaScript constant Math.PI
    var L = new Array();   // Define the array.  It's empty now.  
    for (var i=0; i< Ntrials; i++) {
         var Ntot = (i+1)*Ncoils;  // Total number of coils this sample.
         L[i] = 2*pi*R*Ntot;   // This is the ith value of length in the array
    }

    var resist = [0.6,1.1,1.6,2.1,2.7,3.5];   // Ohms

    // ********************* New code block ***********************************

    var fitObj = doFit(L,resist,1,4); 
    var fitData = fitObj[0];
    var fitEqu = fitObj[1];
    var slope = fitObj[2][0];

    var A = 4.64e-7;    // m^2
    var rho = slope*A;

    fitEqu = 'R = '+fitEqu.replace('t','L');
    // ******************** End new code block ********************************

    // Graphics commands
    optsGraphTitles[0].name = 'Length (m)';    // x-axis title
    optsGraphTitles[1].name = 'Resistance (\\Omega )';  // y-axis title
    optsGraphTitles[2].name = 'Resistance vs Length';   // main title
    optsGraphTitles[3].name = '\\rho = '+rho.toExponential(3)+' \\Omega m; '  // sub title

    // optsLegend = false;   // Old code
    optsLegend = true;      // default
    optsSeries[0].name = 'Data'
    optsSeries[1].itype = 1;      // make it a line
    optsSeries[1].name = fitEqu;  // show the equation of fit  
    xprec = 1;               // Set the precision on the x axis

    // Old code
    //drawSeriesPlot([L],[resist], 0,0, 'resistivityPlot.svg');
    // New code
    drawSeriesPlot([L,L],[resist,fitData], 0,0, 'resistivityPlot2.svg');

The resulting graph looks like this.

 

 

 In the next post, I'll give some more sample codes for plotting fits and expose some more graphics parameters for adjusting your plots' appearances.  Cheers!

My First Script: Plotting Data

This post is the first in a series designed to show you how to use OT/JSvg to analyze data.  There's a whole lot to JSvg, so rather than attempt to cover everything at once, we'll take it piecemeal.  Perhaps the simplest analysis used in many science classes the plotting of data.  OT/JSvg makes this relatively easy to do using regular JavaScript, the linga franca of the internet.  But to get started, you will need to get the Analysis Console open first.  In Firefox, go to the Tools option in the Menu Bar, and select "Start Analysis Console".  If you don't see it, you don't have OT/JSvg properly installed, and I'll refer you to the previous post.   At this point, you will see the console pop up.

Figure 1: A blank console after it has first popped up.

The purpose of the console is to allow users to easily make use of the JSvg libraries and see results in an organized manner.  The console also allows you to open up both the Web Console, a command line interface that allows you to type commands and do calculations using the scope of the window, and ScratchPad, a simple JavaScript programming editor with syntax highlighting.  Today, we're going to use ScratchPad to write the code.  The results will be displayed in the Analysis Console.

Let's start with the data.  In one particular physics lab, students measure resistance of various lengths of wire.   Of course, as wire is conductive, it takes a whole lot of wire to provide a measurable resistance with most multimeters.  So we use coils of wire rather than linear lengths.  Then the students make their measurements of resistance every so many coils and calculate the length based on the number.  If we know the radius of the coil, and the number of coils, it is pretty straight forward to calculate the total length.

    var R = 0.02275;   // m
    var Ncoils = 100;  // Number of coils
    var pi = Math.PI;  // Use pi as a short hand for the JavaScript 
                       // constant Math.PI
    var L = 2*pi*R*Ncoils;   // This is the total length;

But wait! We can do better than this. Let's say that we make six measurements at increments of 100 coils. We could write out each calculation, OR, we could do one better and use a loop and arrays. Loops allow you to repeat calculations an arbitrary number of times without having to type out each one. In short, they are time savers. Arrays allow you to store results. While in general, they are of a related type (for example lengths of wire), in JavaScript arrays can consist of all kinds of unrelated objects. More on that later. Let's change the code to this.

    var R = 0.02275;   // m
    var Ncoils = 100;  // Number of coils
    var Ntrials = 6;   // Number of trials 
                       //   (presume we sample resistance every 100 coils)
    var pi = Math.PI;  // Use pi as a short hand for the JavaScript constant Math.PI
    var L = new Array();   // Define the array.  It's empty now.  
    for (var i=0; i< Ntrials; i++) {
         var Ntot = (i+1)*Ncoils;  // Total number of coils this sample.
         L[i] = 2*pi*R*Ntot;   // This is the ith value of length in the array
    }

The loop is called a for loop and it works like this. Starting with an index (i) of 0, do this loop while i is less than Ntrials (6 times), incrementing i by 1 each time (i++). The result is that the array L will contain six values based on a formula that was executed six times, once per loop. The formula

    Ntot = (i+1)*Ncoils

calculates the total number of coils in that iteration. The first time, since the loop index i is zero, the result will be:

    Ntot = (0+1)*100 = 100 coils.

The second time, when the index is one, the result will be:

    Ntot = (1+1)*100 = 200 coils.

And so on. We will use the L array as our x values. For y values we have our measured resistances.

   resist = [0.6,1.1,1.6,2.1,2.7,3.5];   // Ohms

Note that there are six values in this array and that they are separated by commas. To this point, one could run the ScratchPad code without the help of JSvg since it's pure JavaScript. However to plot these data, we need to use JSvg graphics functions. The simplest scatterplot function available (without going low-level) is

    drawSeriesPlot([2D xarray], [2D yarray], 0,0, filename)

where the first argument is your x data, and the second is your y data. These are 2D arrays rather than 1D arrays by default so that if you have multiple series of data, you can plot them in a single command. The third and fourth arguments of the function are deprecated. For now, just use zeroes every time you call the command. The final argument is the filename in string form. This means that there should be quotes around the name. It is important that you specify the ".svg" extension to ensure that the graph displays properly. So here's the plotting command.

    drawSeriesPlot([L],[resist], 0,0, 'resistivityPlot.svg');

The 1D arrays L and resist are each enclosed in brackets to make them into 2D arrays. If you are starting with 2D arrays, you will not need to do this. When you run this command your graph will appear in the upper left-hand pane of the console (current graph), and when you mouse-over it, its background will turn orange. Clicking on the orange will make it appear in the main frame.

Figure 2: Results thus far

 You will notice a couple of things at this point: first, that there are no labels on the axes, and second, that there is a legend with a series called "unnamed0" displayed.  Let's add the titles via code (and later we'll see how to add them post-coding). The following commands create titles for your graph.

    optsGraphTitles[0].name = 'Length (m)';    // x-axis title
    optsGraphTitles[1].name = 'Resistance (\\Omega )';  // y-axis title
    optsGraphTitles[2].name = 'Resistance vs Length';   // main title
    optsGraphTitles[3].name = 'for a coil of wire';     // sub title

Figure 3: The graph after titles have been added.

If you're sharp, you'll notice that the \\Omega turned into an upper case omega in the y-axis title.  Generally, you can create Greek letters in this manner.  However, you will also need to remember to leave a space after the letter for the parser to work correctly.

We still have the issue of the legend.  One generally doesn't use the legend for single series plots.  So let's remove it.  Set the variable optsLegend to false to suppress its appearance.

    
    optsLegend = false;

The final graph (for today) will look as follows:

Figure 4: The final graph, no legend and with titles.

It is, of course, possible to change things around without dropping to ScratchPad. When your graph is displayed within the Analysis Console, you can click to change most items. Clicking on a title allows you to change its font or position. Clicking on an axis allows you to set the color of the grid lines, the number of divisions on the axis, the use of minor ticks, the precision of the displayed numbers, and so on. The font of the axis matches the font of the axis title. Should you forget to add titles, mouse over the upper left hand side of the graph, and after a couple of seconds, a title dialog will appear. In the meantime, here's the code for the day.

    
    var R = 0.02275;   // m
    var Ncoils = 100;  // Number of coils
    var Ntrials = 6;   // Number of trials 
                   //   (presume we sample resistance every 100 coils)
    var pi = Math.PI;  // Use pi as a short hand for the JavaScript constant Math.PI
    var L = new Array();   // Define the array.  It's empty now.  
    for (var i=0; i< Ntrials; i++) {
         var Ntot = (i+1)*Ncoils;  // Total number of coils this sample.
         L[i] = 2*pi*R*Ntot;   // This is the ith value of length in the array
    }

    var resist = [0.6,1.1,1.6,2.1,2.7,3.5];   // Ohms

    // Graphics commands
    optsGraphTitles[0].name = 'Length (m)';    // x-axis title
    optsGraphTitles[1].name = 'Resistance (\\Omega )';  // y-axis title
    optsGraphTitles[2].name = 'Resistance vs Length';   // main title
    optsGraphTitles[3].name = 'for a coil of wire';     // sub title

    optsLegend = false;

    drawSeriesPlot([L],[resist], 0,0, 'resistivityPlot.svg');

In the next post, we'll talk about how to fit a line to these data. Cheers!

Another Update

One of the hard parts about developing Firefox extensions is that you are developing for a moving target.  The internals of Firefox (which you want to be able to access) are very complex, and often despite trying to keep up with the development roadmap, one still sometimes misses a step.  When you add that to a complicated codebase like OT/JSvg, it's a recipe for trouble every so often.

So last night, I was running my physics lab and everything was going along smoothly.  Students had updated OT/JSvg or installed it for the first time without any issues and were busily collecting data.  Code editing went smoothly: students entered their data into arrays and were plotting the results, when one of the students running Windows noted that the left hand panes where one can click to view recent graphs had stopped working.  After going around to all the machines which I had dutifully updated to FF26.0 and OT/JSvg 0.11.3, I saw that it was true for all the Windows users.  The sole MacOSX person had no such problem nor did I (running Linux). 

The code does OS detection as it loads and sets the variable OS with its findings.  Then later when OS specific code is run, it makes use of this.  It was clear that this was borked.  After some exploring I found that there is new (??) internal object called OS which contains a wealth of info about the operating system.  Therefore, tests as such:

if (OS == 'Windows') {

// Do windowy stuff here

} else {

// Do posix (Linux/MacOSX) stuff here

}

were breaking since OS was an object and not a string.  Having discovered this, I fixed the JSvg part of the code (for a quick release) and am now moving into looking at the OT portion.

While I was at it, I pushed forward some cosmetic changes to the Analysis Console to make things more findable or otherwise legible.  There are now also tooltips for all the buttons.  Here's a preview:

Graph of last night's data with updated GUI.

Anyhow, that's it for now.  Enjoy the 0.11.4 release.

Installing OpenTrack/JSvg

Being a Firefox extension, OpenTrack is relatively easy to install.  If you are already driving Firefox, you can simply point your browser here to download and install the latest edition.  The video below shows the installation process.


 (If you're not seeing a video, here is the link)

Due to Mozilla's very cool extension architecture, once OT is installed, updates are pushed out automatically unless you really don't want this.  Considering that many updates are either bug fixes or new features, it's not clear that you wouldn't want OT to auto update.

A note on this:  If you are like many people and never shut down your browser, you may not get automatic updates.  In which case, you can either click on Tools|Add-ons to get the following tab.

 Near the search bar on the upper right hand side, there is a tool symbol.  This can be clicked to check for updates.  Note that you can also turn automatic updates on or off.

Finally, to start up the Analysis Console, be sure your Firefox menu bar is visible at the top, go into the Tools pulldown menu, and select (in green, except for Mac) "Start Analysis Console."    This will open a window as shown below.

 

The Analysis Console exists to give you access to the statistics and graphics functions embedded in OpenTrack/JSvg.  At this point you have several options:

• If all you want to do is graph something without a lot of fanfare, click on the "Create Plot" button.  
• If you want to do some calculations (that are not persistent), click on the "Toggle Console" button.  The console works as a very nice line calculator. For example, you can type: "a = 2; b=3; c = a*b;" to do a simple math problem.
• If you want to write code, click on the "Open Scratchpad" button.  ScratchPad runs JavaScript code in the context of the window behind it. (The same is true for the Web Console.)  All of the functions and graphics commands of OpenTrack/JSvg are available to you so long as this is true.

 You can also see sample plots and very simple analyses as well as get a sense of some of the functions that are available by pressing the labeled buttons.    In the next post, I will go through some simple examples involving the use of JSvg.

[Edit] One functionality in this release has been curtailed by an internal change in variable.  Mozilla has now turned OS into an object, so my own OS sniffing is sometimes burping as of OT 0.11.3.  I've got a patch in the works, which seems to be working as of now.   I'll post an update tomorrow.

What is JSvg?

In 2005, I was asked by a colleague to write a 2D video analysis package that would work on any operating system and which was freely available to students.  In order to display kinematic results from tracking, a graphing library was needed.  At the time, there was no jQuery, Raphael, svgjs or other Javascript framework for rapid development or graphics.  Nor were there any statistical libraries available.  To meet the requirements of relatively easy development, easy access for users, and since data were to be processed locally, Firefox (XUL+Javascript) was chosen as the platform.  Likewise, the SVG graphical format was chosen since it was possible to create vector graphics that were scriptable and with little work, could be publication quality.  The result of this effort is JSvg.

Fig 1: View of original OT from 2005

JSvg has, to this point, been part of the OpenTrack (OT) Firefox extension.  You won't find OT in Mozilla's addons due to its size (more than 30,000 lines of code), and the need for quick fixes and rereleases mid-semester.  There was no time for the Mozilla approval process.  Rather OT has been released via two sources simultaneously: Sourceforge, and niiler.com (xpi file). The first site was chosen for visibility whereas the second was chosen for automatic updating (which Sourceforge doesn't allow). 

While still embedded as part of OpenTrack, since 2005, JSvg has evolved to become a rather robust graphical and analysis program, run locally as a Firefox extension.  It makes ample use of Firefox's Web Console and Scratchpad and has recently obtained a console for easy navigation of results.

Figure 2: Some of the latest from JSvg showing analysis capabilities comparing different distributions.

So what can it do?  In addition to the many types of plots, it is also capable of running a number of statistics (t-tests, ANOVA, regression, descriptive stats...), signal processing (smoothing, fitting, FFT, derivatives), linear algebra, and more.

The next version of JSvg will be released shortly, both alongside OpenTrack, and for the first time, as its own extension.  This blog will serve as a tutorial to JSvg showing users how to do certain analyses, coding, or graphics.  In a number of cases, I will show benchmarks versus the R statistical language (which is my old standby).  [Note: due to limitations of Javascript, the precision of JSvg results compared to R are not nearly as good.  But that said, the results are often comparable to two or three sig figs.]

Until next time.