Problems with matplotlib._cntr.trace()

[Stephen Terry]

Is anyone else having problems running from source with matplotlib >=
1.0.0? When I load a patient's data, go to the 2D view tab and check
an isodose to display, I'm getting this traceback:

Traceback (most recent call last):
File "C:\Documents and Settings\flounder\Desktop\dicompyler
\baseplugins\2dview
.py", line 426, in OnPaint
self.DrawIsodose(isodose, gc, isodosegen)
File "C:\Documents and Settings\flounder\Desktop\dicompyler
\baseplugins\2dview
.py", line 279, in DrawIsodose
contours = isodosegen.trace(level, points = 0)
TypeError: 'points' is an invalid keyword argument for this function

I poked around in the matplotlib source and it looks like that keyword
might have been removed from the trace function. Here's a link to the
diff

http://matplotlib.svn.sourceforge.net/viewvc/matplotlib/trunk/matplotlib/src/cntr.c?r1=7341&r2=7342

Steve

[Adit Panchal]

Interesting find.

I have not tried it out, because I know matplotlib 1.0.0 has some
funky issues with axes rendering compared to 0.99 and below.
Therefore, I have been using and building the binaries with 0.99. I
wonder if the rendering was resolved in 1.0.1? I will have to take a
look and see.

I could add a conditional so that it calls the function with the
appropriate arguments depending on which version of matplotlib you are
using. I actually just made a conditional like this a few weeks ago to
work around a change in the read_file method in pydicom 0.9.4 vs
0.9.5.

Thanks,

Adit

[flounder.sdt]

Thanks for taking a look at it. I've already tried updating to 1.0.1
and that didn't seem to help.

I'm working on a rudimentary plugin to sum two RTDose objects together
and another to recalculate the DVHs afterward (based on Roy's code from
a few months ago). I'm nearing the point where it would be useful to
get some feedback from other people. What is the best way to distribute
the plugin? Should I just upload the source to this group? Thanks.

Steve

[Adit Panchal]

On Wed, Apr 6, 2011 at 15:55, flounder.sdt <flound...@gmail.com> wrote:
> Thanks for taking a look at it.  I've already tried updating to 1.0.1 and
> that didn't seem to help.

I fixed the issue. Apparently the point argument wasn't even needed in
0.99 and below either as it was optional. The problem was that the
function now returns the vertices and the segments as opposed to
previously only returning the vertices.

Let me know if that works for you. I tested on Mac and Windows with
both 0.99 and 1.0.1.

Unfortunately, I think the axes issue is still present in 1.0.0 and
above. Please see the issue for more details:
http://code.google.com/p/dicompyler/issues/detail?id=39

> I'm working on a rudimentary plugin to sum two RTDose objects together and
> another to recalculate the DVHs afterward (based on Roy's code from a few
> months ago).  I'm nearing the point where it would be useful to get some
> feedback from other people.  What is the best way to distribute the plugin?
>  Should I just upload the source to this group?  Thanks.
>

Regarding your plugin, that is really cool! Hopefully we can get it
integrated somehow.

I had planned on announcing this at some point in the near future, and
probably will make a full announcement in a separate email, but
dicompyler is soon going to have a plugin repository. It will be
located at: http://code.google.com/p/dicompyler-plugins/

Basically, plugin authors can have their code hosted there with any
relevant files or screenshots. Also you can have a wiki page to
document your plugin if you like.

I am still working on the final details on how the repository will be
structured, but what I plan is that each author create a folder with
their plugin with their relevant files, which can include your readme
file. If you prefer, you can create a wiki page for it as well.

After plugins are hosted there, I am planning on updating the plugin
manager in dicompyler to browse the repository and show the available
plugins. Additionally, it will either link to or display the contents
of your readme file / wiki page.

There will be a few changes necessary for the plugin to work on the
repo, but I will post all the relevant info soon.

In the meantime, feel free to upload the source here to the group.

Thanks,

Adit

[flounder.sdt]

As was foretold, here is the source for my dose sum plugin. I've done some
preliminary testing and, so far, it seems to work well. After you open a
patient, you can select "sum RT dose" from the tools menu. This reopens the
Open Patient dialog. If you select a RT Dose object that belongs to the same
image set, the plugin will replace the existing RT Dose in the patient data with
the summed dose.

I would be grateful for any feedback you can give me on the plugin. I see two
main issues right now. The first is accuracy. My testing shows about a 0.3%
error in dose compared to summing the plans in Eclipse. This may be
interpolation error, but I should be able to do better than that.

The second is speed. If the two dose grids coincide, you can just directly sum
the arrays. Otherwise, there has to be some interpolation. If you're running
from source and have SciPy installed, the plugin calls
scipy.ndimage.map_coordinates. Otherwise, it falls back on a trilinear
interpolation function I wrote in Python. My function is incredibly slow.
There is probably some optimization I could do to improve things, but the
easiest thing for me would be to add SciPy as a dependency. I know that's going
to add a lot of bloat to dicompyler, so I understand if you don't want to go
that route.

Anyway, please give my plugin a try and let me know how it works. I'm going to
start getting the DVH calc plugin ready and hopefully I can post that in the
near future.

Steve

plansum.py
-----------

import dicomgui

import wx
from wx.lib.pubsub import Publisher as pub
import numpy as np
import dicom

import os
import unittest

def pluginProperties():
"""Properties of the plugin."""

props = {}
props['name'] = 'Sum RT Dose'
props['description'] = "Adds multiple RT Dose objects together."
props['author'] = 'Stephen Terry'
props['version'] = 0.1
props['plugin_type'] = 'menu'
props['plugin_version'] = 1
props['min_dicom'] = ['rtdose','rtplan']
props['recommended_dicom'] = ['rtdose','rtplan']

return props

class plugin:

def __init__(self, parent):

self.parent = parent

# Set up pubsub
pub.subscribe(self.OnUpdatePatient, 'patient.updated.raw_data')
pub.subscribe(self.OnUpdatePatient, 'patient.updated.parsed_data')

def OnUpdatePatient(self, msg):
"""Update and load the patient data."""
if msg.data.has_key('rtdose'):
self.rtdose = msg.data['rtdose']
if msg.data.has_key('rxdose'):
self.rxdose = msg.data['rxdose']
if msg.data.has_key('structures'):
self.structures = msg.data['structures']
if msg.data.has_key('rtplan'):
self.rtplan = msg.data['rtplan']

def OnDestroy(self, evt):
"""Unbind to all events before the plugin is destroyed."""

pub.unsubscribe(self.OnUpdatePatient)

def pluginMenu(self, evt):
"""Open a new RT Dose object and sum it with the current dose"""
self.ptdata = dicomgui.ImportDicom(self.parent)

if (self.ptdata['rtplan'].SeriesInstanceUID != \
self.rtplan.SeriesInstanceUID):
dlg = wx.MessageDialog(self.parent,
"The image sets for both doses do not match.",
"Cannot sum doses.", wx.OK | wx.ICON_WARNING)
dlg.ShowModal()
dlg.Destroy()
return

old = self.rtdose
new = self.ptdata['rtdose']
sumDicomObj = SumPlan(old, new)
self.ptdata['rtdose'] = sumDicomObj
if self.ptdata.has_key('rxdose') and self.rxdose:
self.ptdata['rxdose'] = self.rxdose + self.ptdata['rxdose']
pub.sendMessage('patient.updated.raw_data', self.ptdata)

def SumPlan(old, new):
""" Given two Dicom RTDose objects, returns a summed RTDose object"""
"""The summed RTDose object will consist of pixels inside the region of
overlap between the two pixel_arrays. The pixel spacing will be the
coarser of the two objects in each direction. The new DoseGridScaling
tag will be the sum of the tags of the two objects."""

#Recycle the new Dicom object to store the summed dose values
sum_dcm = new

#Test if dose grids are coincident. If so, we can directly sum the
#pixel arrays.
if (old.ImagePositionPatient == new.ImagePositionPatient and
old.pixel_array.shape == new.pixel_array.shape and
old.PixelSpacing == new.PixelSpacing and
old.GridFrameOffsetVector == new.GridFrameOffsetVector and False):
print "using direct sum"
sum = old.pixel_array*old.DoseGridScaling + \
new.pixel_array*new.DoseGridScaling

else:
#Compute mapping from xyz (physical) space to ijk (index) space
scale_old = np.array([old.PixelSpacing[0],old.PixelSpacing[1],
old.GridFrameOffsetVector[1]-old.GridFrameOffsetVector[0]])

scale_new = np.array([new.PixelSpacing[0],new.PixelSpacing[1],
new.GridFrameOffsetVector[1]-new.GridFrameOffsetVector[0]])

scale_sum = np.maximum(scale_old,scale_new)

#Find region of overlap
xmin = np.array([old.ImagePositionPatient[0],
new.ImagePositionPatient[0]])
ymin = np.array([old.ImagePositionPatient[1],
new.ImagePositionPatient[1]])
zmin = np.array([old.ImagePositionPatient[2],
new.ImagePositionPatient[2]])
xmax = np.array([old.ImagePositionPatient[0] +
old.PixelSpacing[0]*old.Columns,
new.ImagePositionPatient[0] +
new.PixelSpacing[0]*new.Columns])
ymax = np.array([old.ImagePositionPatient[1] +
old.PixelSpacing[1]*old.Rows,
new.ImagePositionPatient[1] +
new.PixelSpacing[1]*new.Rows])
zmax = np.array([old.ImagePositionPatient[2] +
scale_old[2]*len(old.GridFrameOffsetVector),
new.ImagePositionPatient[2] +
scale_new[2]*len(new.GridFrameOffsetVector)])
x0 = xmin[np.argmin(abs(xmin))]
x1 = xmax[np.argmin(abs(xmax))]
y0 = ymin[np.argmin(abs(ymin))]
y1 = ymax[np.argmin(abs(ymax))]
z0 = zmin[np.argmin(abs(zmin))]
z1 = zmax[np.argmin(abs(zmax))]

sum_ip = np.array([x0,y0,z0])

#Create index grid for the sum array
i,j,k = np.mgrid[0:int((x1-x0)/scale_sum[0]),
0:int((y1-y0)/scale_sum[1]),
0:int((z1-z0)/scale_sum[2])]

x_vals = np.arange(x0,x1,scale_sum[0])
y_vals = np.arange(y0,y1,scale_sum[1])
z_vals = np.arange(z0,z1,scale_sum[2])

#Create a 3 x i x j x k array of xyz coordinates for the interpolation.
sum_xyz_coords = np.array([i*scale_sum[0] + sum_ip[0],
j*scale_sum[1] + sum_ip[1],
k*scale_sum[2] + sum_ip[2]])

#Dicom pixel_array objects seem to have the z axis in the first index
#(zyx). The x and z axes are swapped before interpolation to coincide
#with the xyz ordering of ImagePositionPatient
sum = interpolate_image(np.swapaxes(old.pixel_array,0,2), scale_old,
old.ImagePositionPatient, sum_xyz_coords)*old.DoseGridScaling + \
interpolate_image(np.swapaxes(new.pixel_array,0,2), scale_new,
new.ImagePositionPatient, sum_xyz_coords)*new.DoseGridScaling

#Swap the x and z axes back
sum = np.swapaxes(sum, 0, 2)
sum_dcm.ImagePositionPatient = list(sum_ip)
sum_dcm.Rows = len(y_vals)
sum_dcm.Columns = len(x_vals)
sum_dcm.NumberofFrames = len(z_vals)
sum_dcm.PixelSpacing = [scale_sum[0],scale_sum[1]]
sum_dcm.GridFrameOffsetVector = list(z_vals - sum_ip[2])

sum_scaling = old.DoseGridScaling + new.DoseGridScaling

sum = sum/sum_scaling
sum = np.uint32(sum)

sum_dcm.pixel_array = sum
sum_dcm.PixelData = sum.tostring()
sum_dcm.DoseGridScaling = sum_scaling
return sum_dcm

def interpolate_image(input_array, scale, offset, xyz_coords):
"""Interpolates an array at the xyz coordinates given"""
"""Parameters:
input_array: a 3D numpy array
scale: a list of 3 floats which give the pixel spacing in xyz
offset: the xyz coordinates of the origin of the input_array
xyz_coordinates: the coordinates at which the input is evaluated.

The purpose of this function is to convert the xyz coordinates to
index space, which scipy.ndimage.map_coordinates and trilinear_interp
require. Following the scipy convention, the xyz_coordinates array is
an array of three i x j x k element arrays. The first element contains
the x axis value for each of the i x j x k elements, the second contains
the y axis value, etc."""

use_scipy = True
try:
import scipy.ndimage
except ImportError:
use_scipy = False

indices = np.empty(xyz_coords.shape)
indices[0] = (xyz_coords[0] - offset[0])/scale[0]
indices[1] = (xyz_coords[1] - offset[1])/scale[1]
indices[2] = (xyz_coords[2] - offset[2])/scale[2]
if use_scipy:
return scipy.ndimage.map_coordinates(input_array, indices, order=1)
else:
return trilinear_interp(input_array, indices)

def trilinear_interp(input_array, indices):
"""Evaluate the input_array data at the indices given"""

output = np.empty(indices[0].shape)
x_indices = indices[0]
y_indices = indices[1]
z_indices = indices[2]
for i in np.ndindex(x_indices.shape):
x0 = np.floor(x_indices[i])
y0 = np.floor(y_indices[i])
z0 = np.floor(z_indices[i])
x1 = x0 + 1
y1 = y0 + 1
z1 = z0 + 1
#Check if xyz1 is beyond array boundary:
if x1 == input_array.shape[0]:
x1 = x0
if y1 == input_array.shape[1]:
y1 = y0
if z1 == input_array.shape[2]:
z1 = z0
x = x_indices[i] - x0
y = y_indices[i] - y0
z = z_indices[i] - z0
output[i] = (input_array[x0,y0,z0]*(1-x)*(1-y)*(1-z) +
input_array[x1,y0,z0]*x*(1-y)*(1-z) +
input_array[x0,y1,z0]*(1-x)*y*(1-z) +
input_array[x0,y0,z1]*(1-x)*(1-y)*z +
input_array[x1,y0,z1]*x*(1-y)*z +
input_array[x0,y1,z1]*(1-x)*y*z +
input_array[x1,y1,z0]*x*y*(1-z) +
input_array[x1,y1,z1]*x*y*z)

return output

class PlanSumTest(unittest.TestCase):

def testMethod(self):
ss = dicom.read_file('../testdata/plansum/rtss.dcm')
rtd1 = dicom.read_file('../testdata/plansum/rtdose1.dcm')
rtd2 = dicom.read_file('../testdata/plansum/rtdose2.dcm')

sum = SumPlan(rtd1, rtd2)

self.assertEqual(sum.pixel_array.shape,(4,10,8))
self.assertEqual(sum.PixelSpacing, [5.,5.])

for pos1,pos2 in zip(sum.ImagePositionPatient,[-15.44,-20.44,-7.5]):
self.assertAlmostEqual(pos1,pos2,2)

difference = abs(sum.pixel_array[1,4,3]*sum.DoseGridScaling - 3.006)
delta = 0.01
self.assertTrue(difference < delta,
"difference: %s is not less than %s" % (difference, delta))

difference = abs(sum.pixel_array[1,1,1]*sum.DoseGridScaling - 2.9)
self.assertTrue(difference < delta,
"difference: %s is not less than %s" % (difference, delta))

if __name__ == '__main__':
unittest.main()

[flounder.sdt]

Of course, I post code with a glaring error in it. Line 89 reads

old.GridFrameOffsetVector == new.GridFrameOffsetVector and False):

It should read,

old.GridFrameOffsetVector == new.GridFrameOffsetVector):

A bit of testing that I forgot to take out. Sorry.

Steve

[Anthony Ventura]

If I turn on identical grid summing (remove that "and False") it works
well (for identical grids). Otherwise it never seems to finish on my
PC. I also noticed that if I then anonymize the patient (the only way
I know to save) and re-open it, the Rx dose returns to what it was
before and then the isodose are not normalized correctly.

A great start though, thanks.

Anthony

[Adit Panchal]

I only tried it for identical grids - works great with the
example/test data. However, I have some Pinnacle data that produces
the following traceback:

Traceback (most recent call last):

File "/usr/local/lib/wxPython-unicode-2.8.10.1/lib/python2.6/site-packages/wx-2.8-mac-unicode/wx/_core.py",
line 14614, in <lambda>
lambda event: event.callable(*event.args, **event.kw) )
File "main.py", line 329, in OnUpdatePatientData
self.PopulateIsodoses(patient.has_key('images'), patient['plan'],
patient['dose'])
File "main.py", line 366, in PopulateIsodoses
dosedata = dose.GetDoseData()
File "/Users/apanchal/Projects/python/dicompyler/dicomparser.py",
line 520, in GetDoseData
data['dosemax'] = float(self.ds.pixel_array.max())
File "/Library/Frameworks/Python.framework/Versions/2.6/lib/python2.6/site-packages/pydicom-0.9.5-py2.6.egg/dicom/dataset.py",
line 373, in _get_pixel_array
return self._getPixelArray()
File "/Library/Frameworks/Python.framework/Versions/2.6/lib/python2.6/site-packages/pydicom-0.9.5-py2.6.egg/dicom/dataset.py",
line 368, in _getPixelArray
self._PixelArray = self._PixelDataNumpy()
File "/Library/Frameworks/Python.framework/Versions/2.6/lib/python2.6/site-packages/pydicom-0.9.5-py2.6.egg/dicom/dataset.py",
line 342, in _PixelDataNumpy
arr = arr.reshape(self.NumberofFrames, self.Rows, self.Columns)
ValueError: total size of new array must be unchanged

Haven't gotten a chance to look into it, but if you want I can forward
you the data set if it helps.

One suggestion would be to use threading and dlgProgress from guiutil
with a callback from SumPlan to show the progress of the operation, so
that the interface doesn't stall (similar to the anonymize dialog).

Looking forward to the DVH calc. I had tried to work on it before, but
just haven't gotten the time to get it done. Two things I noticed that
was important for the calc was:

1. Make sure that structure holes are removed from the histogram bins
2. Account for the GridFrameOffsetVector with regards to the
structures as they are usually on the CT planes, while the dose planes
are parallel, but not necessarily the same. I saw that you did use
this in your dose sum code, so I assume you must have already run into
this :)

Regarding SciPy, I have been back and forth on this, but I think at
some point it will be inevitable for it to be included with
dicompyler. It is a great library, but rather heavyweight. If more
plugins start to require it (or the main app can take advantage of
it), we will definitely include it.

Adit

[flounder.sdt]

Thanks for the feedback. The interpolation does take forever if it goes through
the python function. I'm looking at ways of trying to speed that up.

For the second part, I'll have to look into how dicompyler saves the dicom files
and see if I'm not setting the rx dose correctly.

Thanks again!

Steve

[flounder.sdt]

I've gotten that same error before if I wrote to pixel_array in a representation
other than uint32. I put in an explicit conversion to avoid that, so I'm
curious as to what is happening. Do you mind forwarding the data?

>
> One suggestion would be to use threading and dlgProgress from guiutil
> with a callback from SumPlan to show the progress of the operation, so
> that the interface doesn't stall (similar to the anonymize dialog).
>
> Looking forward to the DVH calc. I had tried to work on it before, but
> just haven't gotten the time to get it done. Two things I noticed that
> was important for the calc was:
>
> 1. Make sure that structure holes are removed from the histogram bins
> 2. Account for the GridFrameOffsetVector with regards to the
> structures as they are usually on the CT planes, while the dose planes
> are parallel, but not necessarily the same. I saw that you did use
> this in your dose sum code, so I assume you must have already run into
> this :)
>
> Regarding SciPy, I have been back and forth on this, but I think at
> some point it will be inevitable for it to be included with
> dicompyler. It is a great library, but rather heavyweight. If more
> plugins start to require it (or the main app can take advantage of
> it), we will definitely include it.
>
> Adit

I'll have to educate myself on threading because I've never used it before, but
it does seem like a good solution. Thanks for the suggestion.

Steve

[flounder.sdt]

On 4/7/2011 7:45 PM, Adit Panchal wrote:

> I only tried it for identical grids - works great with the
> example/test data. However, I have some Pinnacle data that produces
> the following traceback:
>

Ok. Your pixel data is uint16, so if I'm going to force it to be uint32 I have
to explicitly set the BitsStored and BitsAllocated tags. Here's a version with
that fixed as well as some other bug fixes.

I've made some small optimizations in the interpolation code, but it's becoming
clear that doing it purely in Python is never going to be fast. I'm exploring
weave now in the hope that I can code it in C.

Steve

plansum.py
----------

import dicomgui

import wx
from wx.lib.pubsub import Publisher as pub
import numpy as np
import dicom

import os
import unittest

def pluginProperties():
"""Properties of the plugin."""

props = {}
props['name'] = 'Sum RT Dose'
props['description'] = "Adds multiple RT Dose objects together."
props['author'] = 'Stephen Terry'

props['version'] = 0.11

return props

class plugin:

def __init__(self, parent):

self.parent = parent

pub.unsubscribe(self.OnUpdatePatient)

old.GridFrameOffsetVector == new.GridFrameOffsetVector):

scale_sum = np.maximum(scale_old,scale_new)

sum_dcm.pixel_array = sum
sum_dcm.BitsAllocated = 32
sum_dcm.BitsStored = 32

sum_dcm.PixelData = sum.tostring()
sum_dcm.DoseGridScaling = sum_scaling
return sum_dcm

def interpolate_image(input_array, scale, offset, xyz_coords):
"""Interpolates an array at the xyz coordinates given"""
"""Parameters:
input_array: a 3D numpy array
scale: a list of 3 floats which give the pixel spacing in xyz
offset: the xyz coordinates of the origin of the input_array
xyz_coordinates: the coordinates at which the input is evaluated.

The purpose of this function is to convert the xyz coordinates to
index space, which scipy.ndimage.map_coordinates and trilinear_interp
require. Following the scipy convention, the xyz_coordinates array is
an array of three i x j x k element arrays. The first element contains
the x axis value for each of the i x j x k elements, the second contains
the y axis value, etc."""

# use_scipy = True
# try:
# import scipy.ndimage
# except ImportError:
use_scipy = False

indices = np.empty(xyz_coords.shape)
indices[0] = (xyz_coords[0] - offset[0])/scale[0]
indices[1] = (xyz_coords[1] - offset[1])/scale[1]
indices[2] = (xyz_coords[2] - offset[2])/scale[2]
if use_scipy:
return scipy.ndimage.map_coordinates(input_array, indices, order=1)
else:
return trilinear_interp(input_array, indices)

def trilinear_interp(input_array, indices):
"""Evaluate the input_array data at the indices given"""

output = np.empty(indices[0].shape)
x_indices = indices[0]
y_indices = indices[1]
z_indices = indices[2]
for i in np.ndindex(x_indices.shape):

x0 = int(x_indices[i])
y0 = int(y_indices[i])
z0 = int(z_indices[i])

x1 = x0 + 1
y1 = y0 + 1
z1 = z0 + 1
#Check if xyz1 is beyond array boundary:
if x1 == input_array.shape[0]:
x1 = x0
if y1 == input_array.shape[1]:
y1 = y0
if z1 == input_array.shape[2]:
z1 = z0

x = x_indices[i]%1
y = y_indices[i]%1
z = z_indices[i]%1

return output

class PlanSumTest(unittest.TestCase):

sum = SumPlan(rtd1, rtd2)

self.assertEqual(sum.pixel_array.shape,(4,10,8))
self.assertEqual(sum.PixelSpacing, [5.,5.])

if __name__ == '__main__':
#unittest.main()
rtd = dicom.read_file('./testdata/rtdose.dcm')
SumPlan(rtd, rtd)

[Adit Panchal]

Great, it seems to work fine now. Do you think you should also set the
High Bit to 31 as well?

According to the DICOM standard:

"Most significant bit for pixel sample data. Each sample shall have
the same high bit.
High Bit (0028,0102) shall be one less than Bits Stored."

Also, maybe we could try using http://pastebin.com for code snippets,
until the repository is ready? You can set the expiration date on the
code to "never", so the link would still be valid in the future. I am
seeing some issues when the code is displayed in Google Groups and in
Gmail (although you can get around it by viewing the original source
of the email).

Thanks,

Adit

[Adit Panchal]

Also forgot, maybe using an attachment instead of code in the body of
the email might be even easier.

Adit

[flounder.sdt]

On 4/8/2011 5:41 PM, Adit Panchal wrote:
> Great, it seems to work fine now. Do you think you should also set the
> High Bit to 31 as well?
>
> According to the DICOM standard:
>
> "Most significant bit for pixel sample data. Each sample shall have
> the same high bit.
> High Bit (0028,0102) shall be one less than Bits Stored."
>

Yes, you're right. Here's another updated version as an attachment.

Steve

[Roy]

I'm travelling right now, but will try this out as soon as I get back.

I'm planning to implement a gamma function plugin for dose map
comparisons, so I'm sure much of your code will be useful for me.

Roy

>  plansum.py
> 14KViewDownload