Category Archives: MATLAB

Visual Acuity and Color Confusion

I’ve been working on a project at work that involves visual acuity. We’re basically asking the question, “How easy is it to notice the difference between two parts of an image with the assistance of color?” As part of the process, we take two grayscale images and make a false-color composite from them. Similar regions of the image are gray (which includes black and white), while the differences show up as one hue or another.

There is an almost infinite number of colors we could use for the false coloring by spinning the color wheel, but it’s easiest to pick a primary color (red, green, or blue) and use its complementary color for the opposite. Some colors are harder to distinguish against bright backgrounds; think about picking out a bright yellow object against a bright white background. While other colors are harder to pick out against dark backgrounds; blue on black, for example. Is there an optimal choice that’s also inexpensive to compute?

Let’s find out. Here’s what happens when you perform the “What’s different?” test using the three obvious choices of red/cyan, green/magenta, and blue/yellow:

(Click any image to see it larger and/or full-size.)

The blue/yellow version is not very good because it suffers from the problems listed above, while the green/magenta and red/cyan false colorings seem to produce fairly equivalent results. What to do? The ever-practical Loren suggested that we see what happens if we ask people with color-blindness color confusion. Turns out, we know a couple of such folks, but they were in meetings. What to do?

Oh, that’s right!! A while back, I wrote some MATLAB code that simulates the two most common forms of anomalous color vision. Let’s run that and see what happens. (The speckling in the images is an artifact of the conversion and doesn’t represent how they actually look to people with protanopia or deuteranopia.)

The results are a no-brainer! (Sorry about the pun.) While the red/cyan and green/magenta images are very similar for the majority of us with normal color vision, the red/cyan images become very difficult to use for people with any kind of anomalous color vision. Even though the green/magenta starts to look a lot like the blue/yellow row, which has its own issues of color discrimination, it makes a very good compromise. Green/magenta it is.

Posted in Color and Vision, MATLAB | Leave a comment

Advocating for the Artificial Pancreas

Editor’s note: A few weeks ago, Kim at Texting My Pancreas let me know that she would be writing a skeptical assessment of the Artificial Pancreas Project (APP), after I expressed some similar skepticism. What appears below is me being the devil’s advocate for the artificial pancreas. While I hold most of these views, I’m not in 100% agreement. Because the APP, if successful, may have an enormous impact on the lives of people with diabetes, Kim and I both think it’s important to weigh the pros and cons seriously. This is just the beginning of the debate. Please chime in by leaving a comment with your thoughts here and on Kim’s opposing viewpoint post.

Conflict of interest statement: I work at a company whose software is used to design, simulate, and (possibly) control various parts of the automated pancreas, but I have no personal involvement in the APP . . . although I’m willing to help out in any way I can. (I also contribute money to JDRF, which is funding the APP.) In addition, among the MathWorks’ various technology offerings are products that can be used for controller design; I make no claims about their fitness for use in a clinical setting. These opinions are my own and may not reflect those of my employer.

Kim has given a very accurate description of the Artificial Pancreas Project (APP) and raised some very valid concerns about it. It’s important to be realistic about what the APP is and is not. We agree that, although it has the ability to change the lives of people with diabetes for the better, it is most definitely not a cure; and it does have a significant amount of risk.

But I don’t think the risks should keep us from moving forward with the APP and seeing it as an important therapy until an actual cure is available. I will lay out why I think it’s important to keep moving forward and trying to get as many people on an AP as safely possible. But first, let’s look at Kim’s three main objections:

  1. She doesn’t want to give up control over her personal treatment actions, especially to a system she doesn’t completely trust.
  2. Neither CGM nor pumps are “fool-proof or absent of errors.”
  3. Faster acting insulins leave little room for error and may lead to more rapid overdose if the algorithms miscalculate or the devices malfunction.

Before I get to each of these, let’s look at why we should welcome an artificial pancreas into our lives. The AP can prevent lows and dampen postprandial excursions. It will get us closer to normoglycemia until researchers develop a true cure to type 1 diabetes. As I understand it, the expectation isn’t perfectly tight “control” — I don’t think that’s realistic to expect from anything other than correctly functioning islets. The goal is to smooth out spikes, eliminate lingering highs, and prevent lows.

Essentially an AP will make the “glucoaster” go away and take most of the burden for balancing insulin and food off the shoulders of people with diabetes (PWDs). We have too many variables to consider when making insulin delivery decisions. In fact, we have heuristics (rules of thumb) instead of rules specifically because we’re treating ourselves without all of the information about the hundred different factors that affect blood glucose (BG), and no one can reasonably expect us to make all of the minute-by-minute decisions required by all those things that affect blood glucose.

Most of us just aren’t as good at achieving our BG targets as we want to be, no matter how earnestly we endeavor to hit them. I’m not being judgmental or suggesting that’s a sign of failure. It’s just a testament to the difficulty of this disease. I think the AP will help reduce the emotional toll of diabetes by helping us get more BG and A1c values that we and our healthcare providers hope to see and by removing a lot of the guilt and despair when our decisions don’t work out as we would like.

It’s natural to be skeptical or hesitant before making big changes, especially when they’re related to our health. But we need to accept our own limitations and fallibility and to accept help from technology when it’s available and reliable.

I don’t want anyone to think that I’m sugar-coating or downplaying the very real risks associated with the AP hardware and with more-or-less removing people from the loop. Every treatment option has risks.* Both multiple daily injection (MDI) and conventional pump therapy have risks. It’s possible to draw up the wrong amount of insulin. Pumps can fail, giving too much insulin or too little. BG meters and CGM sensors can give incorrect information or we can miscalculate the number of carbs in a meal, leading us to incorrect dosing. It’s currently too early to say there is any evidence that the AP has a larger risk exposure than MDI or open-loop pumping.

An AP won’t be any more (or less) fallible than a human pumper. If a pump is going to fail, it’s going to fail. I’ve personally had three pumps fail in about ten years of pumping, one of which led to an overdose of about 40 or 50 extra units of insulin during an infusion set change and two of which led to inoperative pumps. Clearly these are unacceptable events. Every AP pump should be built with even more fault-tolerance and redundancy than the current generation of devices. And it’s likely that some automatic decisions may still require confirmation. For example, if the pump exceeds a certain number of units per hour, the wearer may be required to confirm that they ate or provide an estimate of the number of carbs in the meal.

Based on the kind of questions that the FDA regulators were asking at yesterday’s public workshop, I think the level of testing and validation by the FDA is going to be agonizingly thorough. PWDs will eventually grumble about all the delays getting the AP into our hands. And the legal departments of medical device manufacturers are going to be very reticent to give a thumbs up to any device where “normal use” would cause more risk — and therefore bigger lawsuit costs — than the options currently available.

Furthermore, there’s no requirement that the pumps use the fastest acting insulin available. Patients hopefully will be able to choose whether to use it after consultation with their endocrinologists and gathering feedback from those who are using it. And putting glucose or glucagon in the second chamber of the pump should help reduce the possibility of accidental overdose during normal operation.

Finally, there’s the issue of CGM accuracy. This is the trickiest problem, and the one where I most agree with Kim. I think it’s correct to be skeptical of anything that’s currently in clinical trials — anything that uses CGMS, in fact. Having looked at my own data when I wore a CGM sensor for a week, I saw too many deviations from my BG meter and from how my body felt. Clearly this part of the technology is not sufficient for the fully closed loop that defines the APP.

Does this mean that I’m packing it in? Do I doubt the promise of the whole artificial pancreas project? No. It does mean that I wouldn’t volunteer for a closed-loop trial right now.

But I do think it’s wrong to judge the value of the AP, which is at least 2-3 years away from being submitted to the FDA, using today’s options. What is eventually submitted will have to be much more awesome. Fortunately, CGM technology will continue to improve, and the APP offerings will have to include the most-reliable options or it will likely not be approved.

Is the AP perfect? No. Is it riskier than MDI or open-loop pump therapy. Maybe, maybe not — but that’s not a bad assessment for technology that’s still in the prototype and clinical trials stage. The first regulator-approved offerings (hopefully) will be both more robust and much more accurate.

While it’s important to reiterate that the AP is not a cure — and I sincerely hope that companies selling AP devices never try to market it that way — I think it holds enough promise that I’m willing to set aside my fears about mechanical failures and losing control over the minutiae of a disease that I can’t really control on my own.


* — Remember risk is magnitude of a possible problem combined with the likelihood of that problem occurring.

Posted in 101 in 1001, Diabetes, Health Care, MATLAB, NaBloPoMo, NaBloPoMo 2010 | Leave a comment

An Invitation to Diabetes Researchers

I just finished reading Michael Bliss’s excellent The Discovery of Insulin about the amazing work at the University of Toronto in 1921-1923 by Banting, Best, Collip, and Macloed. For most readers, it’s surely a story about discovery and rivalry and collaboration in medicine, culminating in the first effective treatment of diabetes and (very quickly thereafter) the Nobel Prize for Medicine.

For me it was also a history of what I avoided by being diagnosed with type 1 diabetes after 1922. Bliss includes a small sample of the lives of some of the people afflicted with diabetes before the discovery of insulin. I was truly inspired by those who were fortunate enough to receive this “miracle drug,” but I was heartbroken by all of the people of that age who didn’t make it because the only treatment was to survive on a meager 500-or-so calories for intolerable months until slipping into the coma of ketoacidosis and then (eventually) death.

It’s been almost 89 years since the first successful clinical use of insulin, but we still don’t have a cure. At best, insulin is the key part of a hormone replacement therapy where people with diabetes try to mimic a pancreas. At worst, insulin is a fickle treatment that is difficult to use, expensive, and out of reach of millions of people worldwide.

While I’m so grateful for what I have — a treatment that gives life — what we need now is a cure . . . even if it’s only useful for people who are not yet sick or are newly diagnosed.

If you’re working on a cure or a new therapy for diabetes, the rest of this post is for you.

I will help you with any technical software aspects of your research that you need. If you have questions about MATLAB or C++, I will do my very best to answer them. If you need help with software design or object-oriented analysis/design, I can help you structure your solution. If you require some help with the non-development aspects of software development — such as project management or using tools — I will help you there, too.

My particular areas of expertise are MATLAB, C++, medical imaging formats (such as DICOM), some familiarity with image processing, and multithreaded programming using Intel’s TBB. I also know my fair share about life with diabetes. I can’t promise that I’ll have all of the answers to your particular question, but I’ll work with you as much as I can. And I’ll connect you with my coworkers as much as I can, if it comes down to that. Finding a cure that’s safe, effective, and universally accessible isn’t my day job, but moving into the third era of diabetes where no one dies from it and no one has to impersonate a pancreas is something I’d love to be a small (even anonymous) part of.

You can e-mail me at work or home, and you can hit me up on Twitter, too.

Posted in 101 in 1001, Book Notes, Diabetes, MATLAB, Software Engineering | Leave a comment

Reading a Freestyle Blood Glucose Meter with MATLAB

Here’s some MATLAB code that will read data from a Freestyle blood glucose meter that is attached to your computer via serial cable. (Be sure to replace the “smart” quotes with appropriate single quotes. And the “…”, too. Sorry.)

% Jeff Mather - 3 October 2001 % Connect to meter. s = serial('com1', 'baudrate', 19200, ...            'FlowControl', 'Software'); fopen(s); % Initiate communication. fprintf(s, 'mem'); fprintf(s, 'log250'); % Read diagnostics. fgetl(s);   % First line is empty. meter_ID     = fgetl(s); software_ver = fgetl(s); current_time = fgetl(s); num_readings = str2num(fgetl(s)); for p = 1:num_readings     reading = fgetl(s);     sprintf(reading)         pause(0.25);     end fclose(s);
Posted in Computing, Data-betes, Diabetes, Fodder for Techno-weenies, MATLAB | Leave a comment

Why I Love My Job

So you already know that I ♥ my pancreas. I hope that soon I’ll be able to ♥ both my real and “artificial” pancreases. Kerri at Six Until Me has a great write-up of an announcement between JDRF, Johnson & Johnson, and DexCom which should speed my new love’s arrival. (Update: The Diabetes Mine has even more details and an interview.)

And I ♥ my job. I won’t try to take credit for any of the hard work that JDRF-funded scientists and industry R&D folks are doing to make an artificial pancreas a reality. They deserve all of the credit and so much more. But I know for certain that many people at those institutions and others are using our tools as part of the day-to-day toil of making this a reality. (MATLAB’s graphics have a very *ahem* distinctive appearance, which I’ve seen in presentations about the artificial pancreas and diabetes self-management findings.)

The MathWorks‘ mission statement is that we “accelerate the pace of engineering and science worldwide.” Although it might just sound like nice words on a web site or T-shirt, at times like this it makes me very proud to know that we really do make great things happen faster. And knowing that — even as a peripheral actor — I help with progress, well, that inspires me to do my job that much better, by putting a little more care and polish into my features (especially those related to medical imaging), by rooting out all of the hidden software defects that I might otherwise overlook, and by working hard to get as many of those customer-requested features into the product that I can. That may sound a little corny, but I know that everything I do to accelerate the development of our products helps other people accelerate the creation of their products. And I’m not the only one in the office who feels this way.

Now, back to (indirectly) helping other people fight the good fight. . . .

Posted in Diabetes, MATLAB, Software Engineering, This is who we are | 1 Comment

Diabetes Application – A Sample Report

A couple weeks ago, I presented a teaser of my diabetes application. I’ve been working on it a lot lately, creating GUIs to enter most of the data on my insulin pump and blood glucose meter. But mostly I’ve been creating charts and reports. (Putting data in is much less interesting than analyzing it, in an effort to improve my “control.”)

Here’s one of the reports. It has a bit of everything, including the kitchen sink. The one I’m taking with me to my endocrinologist tomorrow is much more focused. Visualizing all of the data associated with diabetes is kind of tricky, so I’m working on that. More features to come. . . .

(Please, don’t judge me based on my readings. :^) Right now I’m using food to cover my insulin, which is the opposite of how it should work — it seems I never really got switched from the NPH mindset. Getting past that is one of the main reasons why I’m doing in this project.)

I’m committed to releasing everything under some kind of open source license. Everything, that is, except the data and anything that would involve a diagnosis or recommend a particular kind of treatment. I don’t want to be on the hook for FDA approval and all that involves.

Posted in Data-betes, Diabetes, Historical Record, Life Lessons, MATLAB, Software Engineering | Leave a comment

MATLAB and Greenspun’s 10th Rule

Just to get this disclaimer out of the way, let’s remind everyone that this is my personal website. Even though I may write some of the articles at work about things I’ve done at work, they aren’t in anyway supervised or endorsed by my employer. I’m only including this disclaimer because this article is meant to be self-deprecating in that “I had the best of intentions, but now I feel a little silly” way.

Several years ago I read Greenspun’s Tenth Rule of Programming: “Any sufficiently complicated C or Fortran program contains an ad-hoc, informally-specified, bug-ridden, slow implementation of half of Common Lisp.” Basically it’s a joke about feature creep in software. I was reminded of this quotation today when I started reading Joshua Kerievsky’s Refactoring to Patterns. He describes the perils of over-engineering, which is basically the same as wasting money.

You see, about seven or eight years ago, I added a Lisp parser to a feature in MATLAB’s Image Processing Toolbox. If you have the toolbox, just take a look in toolbox/images/iptformats/private/dicom_create_IOD.m. There you’ll see a part of the code that looks like this:

function tf = check_condition(expr, metadata) %CHECK_CONDITION  Determine whether a condition is true. % %   Conditions are LISP-style cell arrays.  The first element %   is a conditional operator, remaining cells are arguments %   to the operator. % %   Conditions can be nested.  Each cell array in expr indicates %   a new conditional expression. %[snip] % % Process conditional expressions recursively. % switch (lower(operator)) case 'and'         % This AND short circuits.     for p = 1:numel(operands)         tf = check_condition(operands{p}, metadata);                 if (~tf)             return         end     end [snip]

So why did I do this? The DICOM file format has conditional elements, which only show up in the file if other elements are present, are absent, or have a particular value. I wanted to create a general-purpose, reusable, extensible mechanism for describing these conditions so that I didn’t actually have to write code for all of those conditionals. Moreover, I wanted our customers to be able to write the conditions in a nice tabular form rather than writing code, because — you know — one day they would want to write more of the dozens of kinds of DICOM information objects from scratch than we supported right out of the box.

You probably see where this is going. We ended up extending the product in a different direction after we learned more about how our customers actually wanted to use the DICOM export functionality. So, if you need a Lisp-like algebraic Boolean logic parser in MATLAB, just go look in the previously mentioned file. (Also take a look at the dicom_modules.m file in the same directory to see how to express those conditions.) But for g-d’s sake, ask yourself if you really, really need that, or if you aren’t just fooling yourself into over-engineering your code.

Posted in Computing, Life Lessons, MATLAB, Software Engineering | Leave a comment

Electronic Imaging: 3D

I am attending Electronic Imaging 2009 in San Jose. Well, today was my last day here. Tomorrow I return to snowy Massachusetts. (It’s currently about 60º here.) But I’ve had my fish taco fix for another year, and I’m excited to get back to Lisa and the kitty.

As usual, I learned something new and talked to a lot of people about imaging and MATLAB. Last year, the new-to-me things concerned image quality for mobile devices and image forensics. This year’s big surprises were 3-D imaging, which apparently is going mainstream. (According to Andrew Woods, the chair of the Stereoscopic Displays and Applications conference, there are already two million televisions in people’s homes which are capable of displaying 3D video.) Every year that I attend, I seem to get a different view of the conference; so I’m cautious to read too much into any given year. But somehow, 3D felt different this time.

The exhibitors behind us were demonstrating their more advanced 3D displays and video capture systems, and the person across from us was selling 3D notecards and prints. So my coworker and I thought we’d get in on the action during the waning hours of the exhibit. We decided to try our hand at making an anaglyph, one of the easier kinds of 3D images to make (along with stereograms, which are a bit lo-fi for MATLAB, after all).

The recipe calls for combining two images of the same scene from slightly different perspectives using a simple function that combines the red channel from the left image and the green and blue channels from the right image. And of course, you need the funny glasses. First, the images that we made with my camera phone:


Left image


Right image

Using some MATLAB code that we found on MATLAB Central, you get the following image:


The anaglyph

(Unfortunately, I can’t help you with the funny glasses. Go get yer own.)

I am told the results are “okay for a cha-cha picture.” I don’t have binocular vision, so what do I know?

Have fun. (And don’t be like the guy who said that the 1.57% of the population with monocular vision aren’t really worth worrying about.)

Posted in Color and Vision, MATLAB, Software Engineering | 1 Comment

Article on High Dynamic Range Imaging and MATLAB

Peter Webb, one of my coworkers, published an article entitled “Rendering High Dynamic Range Images on the Web” in the July 2008 issue of the MATLAB Digest. Enjoy!

Posted in Color and Vision, Fodder for Techno-weenies, MATLAB, Photography | Leave a comment

Quickly Finding Numeric Patterns in MATLAB

Comrade programmers!

About five years ago I wrote a little program to help me find numeric patterns in arrays. Basically I needed to find byte patterns like 0xFFFE 0xE000 in the data from DICOM files. Something like MATLAB’s strfind function, except more general than looking for substrings in strings.

So I wrote a rather naïve algorithm to do just that. It wasn’t very fast, but it worked for exploratory use. But then today I needed to do something like that in production code. So I asked one of my coworkers for tips with a faster implementation.

Turns out the strfind works just as well on numeric arrays as it does on strings. And it’s really fast, too.

% "data" is a 1,200,000 element UINT8 array. % It contains 29 instances of the pattern. tic; offsetTable = strfind(data, [254 255 0 224]); toc % Elapsed time is 0.027390 seconds.

Update: Loren Shure posted an even better solution on her blog.

Posted in Computing, MATLAB | 2 Comments

Using C++ iterators on MATLAB mxArrays

Here’s a little something for the MATLAB lovers out there who also program in C++. The code sample at the end of this post shows how to add a C++ iterator to an mxArray, which is the basic datatype for MATLAB arrays. This makes it harder to walk off the end of an mxArray during linear traversal of all of the elements in an array and (hopefully) improves the readability of code that accesses the data stored inside mxArrays. The code gives a couple simple-minded examples of element-wise traversal (provided the MATLAB array you give it is a UINT8 array).

If you’re familiar with the C++ Standard Library, you can start to use your wrapped mxArray like an STL container. I use the term “like” advisedly, since it’s not actually a container. To reduce memory overhead, I don’t make a copy of the data that gets put into ArrayWrapper. Consequently, it doesn’t model the container concept that element lifetimes are the same as the “container” lifetime. You’re shouldn’t count on being able to use ArrayWrapper with any ole STL function. But then again, you might be able to use it with some; I haven’t tried.

#include "mex.h" #include <cstddef> // Treat an mxArray like a container. // (See Josuttis. "The C++ Standard Library." 1999. p. 219-220) template<class T> class ArrayWrapper {   private:     T      *v;     mwSize thesize;   public:     // Type definitions     typedef T         value_type;     typedef T*        iterator;     typedef const T*  const_iterator;     typedef T&        reference;     typedef const T&  const_reference;     typedef mwSize    size_type;     typedef ptrdiff_t difference_type;     // Constructor     ArrayWrapper(const mxArray *theArray)     {         v = static_cast<T*>(mxGetData(theArray));         thesize = mxGetNumberOfElements(theArray);     }     // Iterator support     iterator begin() { return v; }     const_iterator begin() const { return v; }     iterator end() { return v + thesize; }     const_iterator end() const { return v + thesize; }     // Direct element access     reference operator[](size_type i) { return v[i]; }     const_reference operator[](size_type i) const { return v[i]; }     // Size is constant     size_type size() const { return thesize; }     size_type max_size() const { return thesize; } }; void mexFunction(int nlhs, mxArray *plhs[],                  int nrhs, const mxArray *prhs[]) {     const mxArray *inputArray = prhs[0];     ArrayWrapper<uint8_T> inData(inputArray);     // 1- Count the number of samples in the image.     mwSize count = 0;     ArrayWrapper<uint8_T>::const_iterator i = inData.begin();     while (i++ != inData.end())     {         count++;     }     mexPrintf("The image has %ld samples.\n", count);     // 2 - Halve the image values.     mxArray *outputArray = mxDuplicateArray(prhs[0]);     ArrayWrapper<uint8_T> outData(outputArray);     i = inData.begin();     ArrayWrapper<uint8_T>::iterator j = outData.begin();     while (i != inData.end())     {         *(j++) = *(i++) * 0.5;     }         plhs[0] = outputArray; }
Posted in C, Computing, MATLAB | 2 Comments

MATLAB Performance Tricks #1

Avoid str2double and str2num. Use sscanf instead.

For scalars, you’ll see a modest improvement.

>> str = '0009';
>> tic; for p=1:1000, str2num(str); end; toc
Elapsed time is 0.126388 seconds.
>> tic; for p=1:1000, sscanf(str, '%d'); end; toc
Elapsed time is 0.022299 seconds.
>> str = '3.14159265';
>> tic; for p=1:1000, str2double(str); end; toc
Elapsed time is 0.056466 seconds.
>> tic; for p=1:1000, sscanf(str, '%f'); end; toc
Elapsed time is 0.017805 seconds.

For vectors, you’ll see a more hefty speed up.

>> str = '0009 3.14159265';
>> tic; for p=1:1000, str2double(str); end; toc
Elapsed time is 0.480512 seconds.
>> tic; for p=1:1000, sscanf(str, '%f'); end; toc
Elapsed time is 0.027449 seconds.

Also favor sprintf instead of num2str.

>> tic; for p=1:1000, num2str(p); end; toc
Elapsed time is 0.101599 seconds.
>> tic; for p=1:1000, sprintf('%d', p); end; toc
Elapsed time is 0.018325 seconds.
Posted in Computing, Life Lessons, MATLAB | Leave a comment

Security Vulnerability in CDF, plus a MATLAB Fix

The CDF folks at Goddard Space Flight Center have identified a security vulnerability — a buffer overflow to be specific — that can enable the execution of arbitrary code on your machine if you open a particular malformed file. If you’re accessing CDF files via MATLAB, you can download a security patch from NASA GSFC.

Thank you. That is all.

UPDATE: You can also download an update directly from The MathWorks.

Posted in Computing, File Formats, MATLAB | Leave a comment

Electronic Imaging 2008: Color Universal Design (and a MATLAB-based simulator)

Note: Here’s another dispatches about what happened at the Electronic Imaging symposium even though I’m back from San Jose. They will continue until I run out of useful things to write or something more interesting happens in my life.

Yasuyo Ichihara of Kogakuin University presented some research that she and her colleagues performed on what many of us call “color blindness” or “color deficiency.” She would prefer that we say these differently sighted folks have “color confusion,” which I’ll probably try since it’s slightly more evocative than “color deficiency.” And from a PC point of view, it may not be nice to stagmatize roughly 5% of the population as deficient, etc.

Anyway, they presented research of how their lab applied color universal design to remake the map and timetables of the Tokyo subway system to be friendlier for people with both protanopia and deuteranopia. It turns out that the color confusion can be quite simply modeled by placing two iso-chromaticy lines on the standard x-y chromaticity diagram. Along those lines, all colors appear to be the same hue; so it’s not too hard to pick colors that are maximally contrasting by avoiding the iso-lines. They picked four colors, which they claim are suitable for many print applications: black, orange-red, bluish-green, and a more pure blue.

Color universal design doesn’t stop at picking the right colors. Adding non-color information — such as letters, underlining, framing, or other symbols — proved to be very useful aids to comprehension, too.

All of this reminded me that awhile back I implemented a color deficiency simulator in MATLAB based on work by Hans Brettel. I’m making these functions available here in the hopes that they’re useful and that you actually use them when designing something colorful. Simply download the protanopia and deuteranopia simulator M-files and the associated data file. (You will need the Image Processing Toolbox, since these functions use rgb2ind.)

Let’s see what the peppers image I created some years ago looks like to people with color confusion:

rgb = imread('peppers.png'); figure; imshow(rgb); title('Normal color vision') figure; imshow(deuteranopia(rgb)); title('Deuteranopia') figure; imshow(protanopia(rgb)); title('Protanopia')


Posted in Color and Vision, Computing, MATLAB | 3 Comments

Film-like tone mapping

Every once and a while someone will ask me why I still use film. First of all, I have a film camera, and to buy a new one would be a significant chunk of cash. But that’s not the main reason. I like the look of film; I know how the various emulsions work and know what I need to do to get the best image from them; and (above all) I like the fact that I have an artifact when I get my film back from the lab.

Digital does have a lot of advantages, though. It provides instant feedback at the point of making the photograph. It’s also hard to beat the speed and flexibility of a digital workflow. And for 35mm-sized sensors, it’s hard to tell the difference from film.

This article mixes the look of film and the enormous possibilities of digital. It also continues the trend of notes about high dynamic range images that I’ve been posting recently.

Color film is an analog tone reproduction operator. It renders a vast range of illumination values to a particular color image. That is, it maps the high dynamic range world to a lower dynamic range representation. The major factors that influence this tone mapping “function” include the exposure settings (shutter time, aperture, and ISO speed), the film response curves, and development technique. A film response curve for Fuji Velvia, the color film I use, appears below. Notice how the curve is more-or-less flat in the highlights and shadows and that it has an S-shaped curve in between.

Velvia film response curve

If you take some liberties with that S-shaped part of response curve where the real detail in the image is, you might almost say that the curve has a flat part for the shadows, a flat part for the highlights, and a line that connects the two.

When I first started puttering around with high dynamic range images, I began by translating what I knew about how film renders a scene with a lot of tones:

  • Film is sensitive across a particular exposure range.
  • Parts of the image with less than a certain amount of exposure are rendered as black.
  • Overexposed parts of image are rendered as white (technically clear).
  • You can’t change the film sensitivity (except with some processing tricks that we’ll ignore).
  • You can change the amount of light that is recorded during the exposure.

It wasn’t hard for me to combine this information with the generalization about the film response “curve” to produce a very simple tone-renderer which simulates color slide film in MATLAB. You can download the M-file from my page on MATLAB Central, where you can find some other good stuff. (I’ve also reprinted it below.) I ran it on the HDR image of my office, varying the exposure from 0 EV to 12 EV in half-stop increments, and was very glad to see the same kind of results that I do with film.

My office rendered from HDR to LDR in 25 steps
(Click for larger image…)

But film is so twentieth-century. :-) So I went on to more interesting forms of tone mapping.


function rgbSimulated = simulateFilm(rgbRadiance, nStops) %simulateFilm   Perform film-like tone mapping. %   LDR = simulateFilm(HDR, middleEV) converts the floating-point high %   dynamic range image HDR to a UINT8 low dynamic range image LDR using a %   method that recreates the sensitivity of film.  The middleEV value sets %   the middle exposure value (EV) for the rendering.  Larger EV numbers %   will generate brighter images; smaller values of middleEV will result %   in darker images. % %   The film sensitivity is set to 6 EV, which is comparable to many %   late-model transparency (positive) films. % %   Example %   ------- % %   Render the same HDR image using different "exposure settings." %   Essentially, simulate exposing the same scene from 0 to 12 EV in 1/2 %   stop steps. % %      % Read the HDR image and create a buffer for the rendered images. %      hdr = hdrread('office.hdr'); %      s = size(hdr); %      ldr = ones(s(1), s(2), 3, 25, 'uint8'); % %      % Perform the tone mapping. %      for p = 0.5:0.5:12.5 %          ldr(:,:,:,p*2) = simulateFilm(hdr, p); %      end %      figure; montage(ldr) % Author: Jeff Mather % Copyright 2006-2007 The MathWorks, Inc. % Set the film sensitivity and the mid-tone log-radiance. sensitivity = 6; midPoint = 3; minLogExposure = midPoint - sensitivity/2; maxLogExposure = midPoint + sensitivity/2; % Film works in stops.  Convert radiance to base 2 to compute perception. % Values that are outside the film sensitivity are lost (min or max). rgbRadiance = rgbRadiance * 2^(nStops); rgbSimulated = rgbRadiance; rgbSimulated(rgbRadiance ~= 0) = log2(rgbRadiance(rgbRadiance ~= 0)); rgbSimulated(rgbSimulated < minLogExposure) = minLogExposure; rgbSimulated(rgbSimulated > maxLogExposure) = maxLogExposure; % Convert to RGB. rgbSimulated = uint8(255 * (rgbSimulated - minLogExposure) ./ ...                      sensitivity);
Posted in Color and Vision, Computing, Fodder for Techno-weenies, MATLAB, Photography | Leave a comment