mathjax + gtag

Pages

Tuesday, August 18, 2015

A JTAG/XSVF Library for Arduino

Introduction


I have recently felt the need to incorporate a JTAG port in a project to program a hardware that contained a CPLD. The idea was to both program it and perform some integrity tests on the board. I imagined something using pogo pins, to make it easier and quicker to test everything. I would also write the necessary test routines and generate some kind of report.

With this objective in mind, I have decided to design an Arduino shield to do the job. The testing routines were not really a big deal. And I was sure I would find some JTAG library for Arduino ready to be used. That was not the case.

There were some projects using Arduino to control a JTAG TAP (Test Access Port), but they were all incomplete. And I had no idea what was really JTAG. So I had to study a little bit to make things work for me.

In the end, the challenge proved enlightening. There were some caveats, both from hardware and from software. I'll try to address them in this article.

The library is hosted in here in github. It is also in the new Arduino IDE library manager, so it should be easy to find it.

After all, what is JTAG?


In the internet days, everything starts with a good look at Wikipedia. This is the link to JTAG on Wikipedia. To make a long story short, here is the quote for JTAG definition on Wikipedia:
The Joint Test Action Group (JTAG) is an electronics industry association formed in 1985 for developing a method of verifying designs and testing printed circuit boards after manufacture. In 1990 the Institute of Electrical and Electronics Engineers codified the results of the effort in IEEE Standard 1149.1-1990, entitled Standard Test Access Port and Boundary-Scan Architecture.
What it means is that the name JTAG originally meant an electronics industry association. But after IEEE published the Standard 1149.1, which referred to a "Standard Test Access Port", this test access port (TAP) has become more or less a synonym for JTAG.

The physical layer consists of 5 signals:

  • \(\hbox{TCK}\) - Test Clock input
  • \(\hbox{TMS}\) - Test Mode Select
  • \(\hbox{TDI}\) - Test Data Input
  • \(\hbox{TDO}\) - Test Data Output
  • \(\overline{\hbox{TRST}}\) - Test Reset Input
The \(\overline{\hbox{TRST}}\) signal is optional.

JTAG enabled devices can be all connected together. The signals TCK and TMS (and \(\overline{\mbox{TRST}}\), if present) should be connected to all devices. As a consequence, all state machines in all devices will be always at the same state. The TDI of the JTAG interface is connected to the TDI of the first device. The TDO of a device should be connected to the TDI of the next device. Finally, the TDO of the last device is connected to the JTAG interface. This way, data can be shifted in or out of a big shift register, that gets bigger if you add more devices.


The IEEE standard defines not only what the electric cable signals are, but also defines the logic state machine that must understand these signals. Objectives were to keep the number of signals on the cable as low as possible, but keeping the architecture versatile enough to perform any task.

The standard goes on to define the basics of the Boundary-Scan Architecture, which seems to be what they had primarily in mind at that time. It is a way to of testing your hardware "on the fly", i.e., while the circuit is operating. Boundary-Scan compatible devices make it possible to control the logical value of output and input pins of your integrated circuits to see how the whole hardware will respond to those stimulus.

What is JTAG used for?


It didn't take much to start "abusing" JTAG to perform other things besides boundary-scan. For example:

  • Customized hardware testing for quality control
  • CPLDs and FPGAs programming
  • Microcontroller debugging
In order to understand how these tasks can be acheived, we must understand how the JTAG TAP works.

The TAP


The TAP is a synchronous finite state machine. The following diagram shows the state diagram of a JTAG TAP. Transitions are controlled by the state of TMS on the rising edge of TCK.


The state machine is actually simple. In what follows, xR means DR or IR.
  • Actions on the test logic can happen either on the risign or on the falling edge of TCK.
  • At any time, there are two registers that you have to be concerned with: the Instruction Register (IR) and the Data Register (DR). The actual content of DR depends on the value loaded on IR. IR must have at least two bits.
  • There is a reset state that is easily reacheable. Starting from any state in this diagram, if you hold TMS high for 5 consecutive clock rising edges, the TAP is guaranteed to enter TEST-LOGIC-RESET. When this state is entered, the IR is loaded with either the IDCODE instruction or the BYPASS instruction, so that when the TAP moves into RUN-TEST/IDLE no action will occur. If the signal TRST is present, it assynchronously forces the TAP into TEST-LOGIC-RESET.
  • RUN-TEST/IDLE is the sate where you will usually wait or pass between operations. It is one of the stable states, meaning if you hold TMS on a determined value, you stay in that state. For example, the instruction RUBINST causes a self-test of the system on this state.
  • The two vertical state columns perform similar funtions. The first is meant to read and write to the DR while the second does the same to the IR. The DR is actually a group of registers. The one you actually access is dependent on the contents of the IR, so you can think of IR as a selector.
  • SELECT-xR-SCAN: Temporary state to enter xR operations or just to pass on to the next state.
  • CAPTURE-xR: The current selected DR is loaded on the shift registers. In the case of the IR, this state loads the fixed binary pattern "01" on the bits closer to TDO. The other bits may load other design-specific data.
  • SHIFT-xR: In this state, the data loaded in the shift-register is both serially shifted in xR through TDI and shifted out through TDO.
  • EXIT1-xR: This is temporary state that provides a way to bypass PAUSE-xR and EXIT2-xR and go straight to UPDATE-xR.
  • PAUSE-xR: Nothing happens, the controller is paused while in this state.
  • EXIT2-xR: A temporary state that provides a way to return to SHIFT-xR.
  • UPDATE-xR: Data is latched into xR on the falling edge of TCK.
In other words, we write to IR to tell the device what we want, and write or read to DR to set a property or get a response.

BYPASS is an instruction that turns DR into a one bit register that always capture a zero.

There is no such thing as shifting-in without shifting-out or vice-versa. We always do both.

Playing around with JTAG


The Inpact software from Xilinx has a nice interface that allows you to play with the TAP:


Select "Boundary-Scan", then "Initialize Chain" and then on the Debug menu choose "Enable/Disable Debug Chain". It is rudimentary and not very practical, but has its uses. Hand collecting of the bits shifted in and out of IR and DR can be tedious. That is why a second interface is provided that is slightly more usefull.

For example, on the previous screen, click on "Test Logic Reset", and this state becomes green. If the IEEE standard is folowed by this device, then the "IDCODE" instruction (if present) must have been loaded into IR. If we shift out DR, then we must be able to access this IDCODE. To do that, fill the box next to SCAN DR with 32 zeroes, click on "Execute" and lets see what comes out:
TDO Capture Data: 00000110111001011110000010010011
If we break this binary string into hexadecimal, we get something more meaningfull. This is 06E5E093, which is the IDCODE for XC2C64A.

At this point, I should mention a very usefull file to have in hand. It is the BSDL file for the device. BSDL stands for "Boundary-Scan Definition Language", its a subset of VHDL that defines the boundary-scan parameters of a device. BSDL has been defined on the same IEEE Std 1149.1. There we can look at the instructions that our device will accept. The relevant part for us now is:

... 
attribute INSTRUCTION_LENGTH of xc2c64a : entity is 8;
attribute INSTRUCTION_OPCODE of xc2c64a : entity is
      "INTEST (00000010)," &
      "BYPASS (11111111)," &
      "SAMPLE (00000011)," &
      "EXTEST (00000000)," &
      "IDCODE (00000001)," &
      "USERCODE (11111101)," &
      "HIGHZ (11111100)," &
      "ISC_ENABLE_CLAMP (11101001)," &
      "ISC_ENABLEOTF (11100100)," &
      "ISC_ENABLE (11101000)," &
      "ISC_SRAM_READ (11100111)," &
      "ISC_SRAM_WRITE (11100110)," &
      "ISC_ERASE (11101101)," &
      "ISC_PROGRAM (11101010)," &
      "ISC_READ (11101110)," &
      "ISC_INIT (11110000)," &
      "ISC_DISABLE (11000000)," &
      "TEST_ENABLE (00010001)," &
      "BULKPROG (00010010)," &
      "ERASE_ALL (00010100)," &
      "MVERIFY (00010011)," &
      "TEST_DISABLE (00010101)," &
      -- "STCTEST (00010110)," &
      "ISC_NOOP (11100000)";
 
attribute INSTRUCTION_CAPTURE of xc2c64a : entity is "XXXXXX01" ; 
attribute IDCODE_REGISTER of xc2c64a : entity is "XXXX0110111001011XXX000010010011";
attribute USERCODE_REGISTER of xc2c64a : entity is "XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX"; 
...
 That tells us that the IR has 8 bits and what codes it does understand. Lets try IDCODE again, now the hard way. Put "0000001" in the "Scan IR" box, then click "Execute". Then, if it is not already there, put "00000000000000000000000000000000" in the "Scan DR" box and click the second "Execute" button. You should get
TDO Capture Data: 00000101
TDO Capture Data: 00000110111001011110000010010011
The first line with 8 bits seem to be the mandatory "01" that IR will always load, plus six bits that have a meaning set by the manufacturer. This is what is implied by the line 'attribute INSTRUCTION_CAPTURE of xc2c64a : entity is "XXXXXX01" ;'  Notice that since we must enter "CAPTURE-IR" before "UPDATE-IR", there is no way for the device to guess what instruction is going to be loaded. So this value is either a constant or dependent upon an internal state of the device. We will see later that it has something to do with checking for read/write protect.

The second line with 32 bits is exactly what we got before when we performed the "Shift-DR" straight after "TEST-LOGIC-RESET".

It is possible to use Impact to check that we have really got the right id code, if you go out of "Debug Mode", click on the device and then on "Available Operations" click "Get Device ID":
INFO:iMPACT - Current time: 8/14/15 5:01 PM
Maximum TCK operating frequency for this device chain: 33000000.
Validating chain...
Boundary-scan chain validated successfully.
'1': IDCODE is '00000110111001011110000010010011'
'1': IDCODE is '06e5e093' (in hex).
'1': : Manufacturer's ID = Xilinx xc2c64a, Version : 0
It is possible to crak a little bit more the of the id code. Bits 31 to 28 "0000" are the version number, bits 27 to 22 (6E5E in hexa) are the part number and the final 11 bits "093" (in hexa) are the manufacturer's id and in this case mean Xilinx. The first bit, bit 0 is in fact required to be '1' if an IDCODE instruction is present. Remember that upon TEST-LOGIC-RESET, IR gets loaded with IDCODE or BYPASS. But BYPASS is required to load a '0' at the start of the scan cycle, so this bit is used to identify which instruction has been used upon reset.

Playing around with SVF


Now we know enough about JTAG, lets see what else we can do. Suppose you want to describe a set of operations to be performed on a JTAG TAP controller, for example, suppose you want to program a certain device. How do you describe what has to be done? The answer is a programming language. SVF stands for "Serial Vector Format", and is a file format that specifies how and which boundary-scan vectors should be transferred to a device, and also which should be some of the expected results, like ID code or checksum.

The problem with SVF is that it is too verbose. Good for humans to read, but excessive for computers to deal with. So Xilinx has come up with XSVF, which is a binary form o SVF.

Fortunately for us, Impact is able to generate both SVF and XSVF. In our previous example, go to menu "Output -> SVF File -> Create SVF File..." and choose a name for the file to save data in. Double click on "Get Device ID" then go to "Output -> SVF File -> Stop Writing to SVF File".

You should get a file like this:

There is a quick reference for SVF and XSVF in the apendices of the document XAPP503 - SVF and XSVF File Formats for Xilinx Devices (http://www.xilinx.com/support/documentation/application_notes/xapp503.pdf).

The most relevant information is in the instructions SIR and SDR. These are a version of the Impact gui interface that we have previously used. But besides telling the SVF machine what it should shift in, it also tell it what it should get shifted out. It also specifies masks, so that the relevant bits may be checked. For example:
SIR 8 TDI (01) SMASK (ff) ;
SDR 32 TDI (00000000) SMASK (ffffffff) TDO (f6e5f093) MASK (0fff8fff) ;
This means "shift the following 8 bits into IR: 00000001, all of which are relevant, then shift 32 zeros into DR, all of which are relevant, then compare the value received with f6e5f093, after masking both with 0fff8fff". Sounds familiar? Pretty much what we did by hand before, except the masking part.

The comments in the file are interesting. And we can see a lot of redundant operations. The ID code is checked three times, and the read/write protect is checked twice. In the end, the device is put on BYPASS.

Playing around with XSVF


If we use Impact to generate a XSVF version of the ID code check program, we get this:


A completely binary file, improper for humans to read but good enough for computers. There you have the exact same program as before, but coded in binary XSVF instructions, rather than SVF. I wrote a disassembler in python for debugging, and the output is the following:


We can see the same redundancies that were present in the SVF file.

The previously mentioned XAPP503 has the documentation for the XSVF instructions. Which boils down to shifting stuff into IR or DR, reading stuff back and pulsing TCK while in some state, usually in . Not new stuff, after all we have been through.

Arduino, at last


Now that we know what JTAG is, we can program a microprocessor to do the job for us. Arduino is a nice choice for a number of reasons I don't need to get into. I have searched for a library that would just do what I wanted, but none seemed to work with my hardware. There were several problems, I'll try to address them.

The first issue I found was memory usage. It is absolutely essential that you keep your strings within flash. Even more when you are in debug mode, where you want lots of output to understand what is going on.

The second problem was the Arduino serial interface, which has no flow control. That means that we must provide one. What I did was to always use a fixed block transfer size. I have used the Arduino software serial buffer avoiding unecessary copying to spare memory. There are 64 bytes in this buffer, but only 63 are useable because of the circular buffer implementation. I have managed to increase its size to 256 bytes, the process is documented in the code. I have spent quite some time to make the transfer as fast as possible.

The cable was another issue. Don't make it too long and do some termination/impedance matching. Reflections in TCK will kill any attempt to program a JTAG device. I have used three voltage dividers to convert from the 5 Volts logic to the 3.3 Volts logic of the devices I was using, and I have chosen the resistor values to kinda match the expected impedance without killing Arduino's ATMEGA 328p output drivers.

The best solution for the signal conversion would have been to use buffers. Usually the JTAG cable will have a VCC signal that can be used to sense the device's operating voltage. That does not exclude the impedance matching, that in this case could be done with some series resistors.

The VCC signal of the JTAG cable is used in this project to detect that the cable is actually connected to some hardware.

Here is a picture of the two programmers I have used, the Xilinx one and the Arduino. Each of them is connected to a XC2C64A breakout board from Dangerous Prototypes.


Here it is a copy/paste of the terminal screen while programming an example:

$ ./xsvf ../xsvf/XC2C64A/VHDL-CPLDIntro3LEDinverse.xsvf
File: /home/user/sketchbook/arduino/libraries/JTAG/extras/xsvf/XC2C64A/VHDL-CPLDIntro3LEDinverse.xsvf
Ready to send 22846 bytes.
IMPORTANT: Free memory: 771 bytes.
Sent:    22846 bytes,        0 remaining ()
IMPORTANT: XCOMPLETE
IMPORTANT: ********
IMPORTANT: Success!
IMPORTANT: ********
IMPORTANT: Processed 1417 instructions.
IMPORTANT: Checksum:  0x36/22846.
IMPORTANT: Sum: 0x0033D4CA/22846.
Received device quit: Exiting!
  Expected checksum:  0x36/22846.
  Expected sum: 0x0033D4CA/22846.
Elapsed time: 4.31 seconds.

Now what?

Besides programming devices, it is now possible to use JTAG to actually communicate with devices. One possibility is to create a JTAG TAP in VHDL and use it to control your device.

Another possibility is to use the Arduino JTAG to hack into hardware. It is not unusual to find "lost" and undocumented JTAG interfaces on several devices, if you search a little around the internet, you will have some ideas of what you can do.

Hope you like this, happy jtagging!

Editions


  1. There is a nice tutorial on YouTube int the channel "Electronica y Linux": https://www.youtube.com/watch?v=x5eW9H1GoeA



Sunday, October 19, 2014

How to post code on Blogger using Gist

The problem

How to post nice pieces of code to a blog? The desired properties are:
  1. Beauty.
  2. Language awareness.
  3. Easy to copy to the clipboard.
  4. Easy to integrate in Blogger.
  5. Easy to maintain

One possible solution

Use Gist.You must have an account on GitHub. Gists are versioned text files, so you can alter them and fix any problems on the Gist page without having to deal with your blog page.

How to integrate it into Blogger

After creating your Gist, copy the "Embed URL" to the clipboard:
The link is something like this:
<script src="https://gist.github.com/mrjimenez/bceea54fdf66a9aec7b2.js"></script>
Reserve a space to integrate the Gist in your blog post. Then choose "HTML" in the Compose/HTML radio button. Insert that code in the proper place in your page, and your're done:


Take a look at the snapshot:


Saturday, October 18, 2014

Uploading photos to Facebook using google-chrome on linux

The problem

Uploading photos to Facebook seems to be unresponsive, I get to this screen and when I press the "Choose File" button, nothing happens:

The solution

Go to a file manager (e.g, konqueror) drag the file with the mouse and drop it over one of the buttons. The file will then be accepted:

If you want to unselect a photo, click the corresponding "Choose File" button. Apparently, its intended function is not to choose a file, it is to reset the file to "No file choosen":

After that, click "Upload Photos" and wait, because it takes some time.

Thursday, October 16, 2014

Cuda 6.5 on OpenSuSE 12.3

NVidia Drivers

Make sure you have the official NVidia drivers installed in your system:
  • Run Yast.
  • Click in "Software Repositories"
  • Click in "Add"
  • Choose "Specify URL", then "Next"
  • Repository Name: "nvidia", URL: http://download.nvidia.com/opensuse/12.3/
  • Confirm
  • Go back to Yast
  • Click in "Software Management"
  • Search for "nvidia"
  • Add the following packages: 
  1. x11-video-nvidiaG03-340.46-30.1.x86_64
  2. nvidia-gfxG03-kmp-desktop-340.46_k3.7.10_1.1-30.1.x86_64
  3. nvidia-computeG03-340.46-30.1.x86_64
  4. nvidia-settings-325.15-1.3.x86_64
  5. nvidia-glG03-340.46-30.1.x86_64
  6. nvidia-texture-tools-2.0.6-36.2.x86_64
  7. nvidia-uvm-gfxG03-kmp-desktop-340.46_k3.7.10_1.1-30.1.x86_64
Notice that the above assumes your board is supported by the G03 kernel driver and that you are using the "kernel-desktop". Make sure you choose the proper driver for your board and the kernel driver corresponding to your kernel.

CUDA Installation

Install the CUDA repository. Although the repository is for OpenSuSE 13.1, it will work perfectly with 12.3.
  • Click in "Add"
  • Choose "Specify URL", then "Next"
  • Repository Name: "cuda", URL: http://developer.download.nvidia.com/compute/cuda/repos/opensuse131/x86_64
  • Confirm
  • Go back to Yast
  • Click in "Software Management"
  • Search for "cuda"
  • Add the following packages (some of them will be automatically added): 
  1. cuda-documentation-6-5-6.5-14.x86_64
  2. cuda-cudart-6-5-6.5-14.x86_64
  3. cuda-cufft-dev-6-5-6.5-14.x86_64
  4. cuda-repo-opensuse131-6.5-14.x86_64
  5. cuda-visual-tools-6-5-6.5-14.x86_64
  6. cuda-cufft-6-5-6.5-14.x86_64
  7. cuda-npp-dev-6-5-6.5-14.x86_64
  8. cuda-curand-dev-6-5-6.5-14.x86_64
  9. cuda-license-6-5-6.5-14.x86_64
  10. cuda-runtime-6-5-6.5-14.x86_64
  11. cuda-misc-headers-6-5-6.5-14.x86_64
  12. cuda-samples-6-5-6.5-14.x86_64
  13. cuda-curand-6-5-6.5-14.x86_64
  14. cuda-toolkit-6-5-6.5-14.x86_64
  15. cuda-cublas-6-5-6.5-14.x86_64
  16. cuda-cusparse-dev-6-5-6.5-14.x86_64
  17. cuda-drivers-340.29-0.x86_64
  18. cuda-cudart-dev-6-5-6.5-14.x86_64
  19. cuda-npp-6-5-6.5-14.x86_64
  20. cuda-command-line-tools-6-5-6.5-14.x86_64
  21. cuda-cusparse-6-5-6.5-14.x86_64
  22. cuda-6.5-14.x86_64
  23. cuda-core-6-5-6.5-14.x86_64
  24. cuda-cublas-dev-6-5-6.5-14.x86_64
  25. cuda-driver-dev-6-5-6.5-14.x86_64
  26. cuda-6-5-6.5-14.x86_64
  • Click "Accept".

Testing





Usefull Links

Monday, October 6, 2014

Moving Averages Using Google Apps Scripts

Abstract

A moving average (MA) process, also known as weighted moving average (WMA) is a type of signal filtering that consists in performing a weighted average over a finite sequence of past samples of the original signal. Implementing such scheme in a worksheet is not always straightforward due to the handling of missing values. This article proposes an implementation of such filters on Google Sheets using Google App Scripts, a Javascript like language. Exponential Moving Average (EMA) filters are an important special case of WMA, so they have been implemented on top of WMA.

Keywords

Moving Average, MA, Weighted Moving Average, WMA, Exponential Moving Average, EMA, Noise Removal, Filtering, Google Sheets, Spreadsheets, Google App Scripts, Javascript.

Introduction

In previous articles [1][2], implementations of WMA and EMA have been proposed using the normal spreadsheet function infrastructure. These implementations try to be compatible with existing spreadsheet standard functions so that they can be easily ported to other spreadsheet e.g., OpenOffice.

The inconvenient in this solution is that the formulas are large and as a consequence, hard to read and maintain, making it easy to slip subtle errors. Also, the formula must be called once for each line to be calculated, and this process has a big overhead. The ideal solution should make a single call to a function that would return an array of processed data. The following implementation addresses both issues.

Implementation

Results

The resulting spreadsheet shows a comparison of the previously posted methods with this scripted based one.

Conclusion

Two javascript functions have been developed to implement the missing data weighted moving averages of previous articles [1][2]. The results have been shown to be identical. The javascript based method has the advantage of being much cleaner to maintain, typically requiring a single cell on the spreadsheet.


Wednesday, September 24, 2014

Weighted Moving Averages on Google Sheets

Abstract

A moving average (MA) process, also called weighted moving average (WMA) is a type of signal filtering that consists in performing a weighted average over a finite past sequence of samples of the original signal. Implementing such scheme in a worksheet is not always straightforward due to the handling of missing values. This article proposes a convenient way to implement such filters on Google Sheets.

Keywords

Weighted Moving Average, MA, WMA, Noise Removal, Filtering, Google Sheets, Spreadsheets.

Introduction

The moving average signal is defined by:
\[y(j) = \sum_{i=0}^n w(i) \cdot x(j-n+i)\]
The spreadsheet function \(AVERAGE()\) performs this computation in the particular case where all \(w(i)\) are equal to one. In this case, it is very easy to deal with missing values, as the only thing you need to do is count the number of non-missing values and divide the sum by this count.

When the \(w(i)\)'s are different numbers, it is not enough to count the number of non-missing values, one needs to perform the final normalization using the sum of the correspondent weight coefficients.

Proposed Solution

The proposed solution is to use the function \(SUMPRODUCT()\) twice, first to get the weighted sum, and second to get the sum of the coefficients that have multiplied non-missing data.

Assume that:
  • The spreadsheet has a page called "Filter1".
  • The averaging coefficients are on column Filter1!A.
  • Cell Filter1!B1 has the formula "COUNT(A:A)", which will count the number of averaging coefficients.
  • The current cell is B22.
Then the claim is that the following formula will calculate the correct value of WMA:
\[\begin{array}{l}
=\\
SUMPRODUCT( \\
\qquad OFFSET(B22, -Filter1!$B$1 + 1, 0, Filter1!$B$1, 1), \\
\qquad OFFSET(Filter1!$A$1, 0, 0, Filter1!$B$1, 1)) \\
/ \\
SUMPRODUCT(\\
\qquad ARRAYFORMULA(\\
\qquad \qquad N(ISNUMBER(\\
\qquad \qquad \qquad OFFSET(B22, -Filter1!$B$1 + 1, 0, Filter1!$B$1, 1)))), \\
\qquad OFFSET(Filter1!$A$1, 0, 0, Filter1!$B$1, 1))
\end {array}\]
The details of this expression are as follows:
  • \(Filter1!$B$1\) is \(n\).
  • \(OFFSET(current\_cell,−n+1,0,n,1)\) is used to produce a range of n cells, of which the current cell is the last one.
  • \(OFFSET(Filter1!$A$1, 0, 0, Filter1!$B$1, 1)\) are the weighting coefficients.
  • \(ARRAYFORMULA(N(ISNUMBER(OFFSET(\cdots))))\) applies the function \(N()\) to the boolean result of \(ISNUMBER()\) for each cell in the current range, which will produce an array of zeroes where there is missing data and ones where there is data.

Results

The formula has been tested against some weighting data. The resulting spreadsheet has a plot of the original data, along with the \(AVERAGE()\) data and WMA data for comparison.

Conclusion

A spreadsheet formula for the correct calculation of a weighted moving average has been derived and successfully tested on Google Sheets. The proposed formula deals with missing values in a way similar to the \(AVERAGE()\) function, avoiding the distortions that would be caused either by using zero in place of the missing values or by packing the original series.

Tuesday, September 23, 2014

Exponential Moving Averages on Google Sheets

Abstract

Exponential moving averages (EMA) is a way to remove noise from data series. Unfortunately, useful straightforward spreadsheet support for EMA is absent. This article examines the problems involved and proposes a one line formula solution to add EMA to a spreadsheet.

Keywords

Exponential Moving Average, EMA, Noise Removal, Filtering, Google Sheets, Spreadsheets.

Introduction

Traditionally, implementing an exponential moving average has been done in spreadsheets using a recursive formula, i.e., an auto-regressive (AR) process or an infinite impulse response (IIR) filter. The formula is the following:
\[y(n) = \alpha \cdot x(n) + (1 - \alpha) \cdot y(n-1)\]
There are some problems with this approach:
  1. Strictly speaking, this is not a moving average. Moving average assumes a finite length sliding window under which data is being weighted. There is a growing window, not a sliding window.
  2. The last item implies that the average takes into account every single data sample. The previous formula actually implements an IIR filter. This kind of averaging never forgets a value, although old values certainly get irrelevant over time. It would be nice to have control over the window where the exponential average takes place, especially when one intends to use a small window.
  3. It does not deal properly with missing values. In a spreadsheet, it is common to have missing data values and the previous approach does not allow them. Converting missing values to zero would unacceptably distort the average value. If one insists in applying this formula to the series by packing the original series, the weights applied to the values in the averaging process will not reflect the actual distance in time that these samples might have.
The function AVERAGE(range) is able to deal with this problem quite simply because it uses the same average on every sample, so it is just a matter of dividing the sum by the number of non-blank entries. In a non-uniform averaging like EMA, we need to keep track of which weights were really applied or not due to missing values, and fix the normalizing factor accordingly.

Proposed solution

The function SERIESSUM(a, n, m, x) is defined as
\[SERIESSUM(a, n, m, x) = \sum_{i=0}^n x_i a^{n+m i}.\]
The proposed solution is to use this function twice, first to calculate the weighted sum, and then a second time to calculate the sum of the weights where data is not missing.
Assume that:
  • Cell F1 contains the geometric progression ratio \((\alpha)\);
  • Cell F2 contains the window size \((n)\);
  • Column B contains the raw data;
  • The current cell is B22.
Then the claim is that the following formula will calculate the correct value of EMA:
\[\begin{array}{l}
=\\
SERIESSUM($F$1, $F$2, -1, \\
\qquad ARRAYFORMULA(N(OFFSET(B22, -$F$2 + 1, 0, $F$2, 1))))\\
/ \\
SERIESSUM($F$1,  $F$2, -1, \\
\qquad ARRAYFORMULA(\\
\qquad \qquad N(ISNUMBER(OFFSET(B23, -$F$2 + 1, 0, $F$2, 1)))))
\end{array}\]
The details of this expression are as follows:
  • \(OFFSET(current\_cell, -n+1, 0, n, 1)\) is used to produce a range of \(n\) cells, of which the current cell is the last one.
  • \(ARRAYFORMULA(N(OFFSET(\dotsc)))\) will apply the \(N()\) function to each element of the argument range to generate a new range with zero values in the missing data cells. Without this trick, \(SERIESSUM()\) would use non-missing values as if they were contiguous.
  • \(ARRAYFORMULA(N(ISNUMBER(OFFSET(\dotsc))))\) will generate a range composed of ones where data is not missing and zeros where data is missing.

Results

The formula has been tested against some weighting data. The resulting spreadsheet has a plot of the original data, along with the \(AVERAGE()\) data and EMA data for comparison.

Conclusion

A spreadsheet formula for the correct calculation of an exponential moving average has been derived and successfully tested on Google Sheets. The proposed formula deals with missing values in a way similar to the \(AVERAGE()\) function, avoiding the distortions that would be caused either by using zero in place of the missing values or by packing the original series.