medical image analysis SVN (old)

#!/usr/bin/env python
#
# Copyright (c) Madrid 2008 
# BIT, ETSI Telecomunicacion, UPM
#
# This program is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation; either version 2 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, write to the Free Software
# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
#
#
#


modules = {'Section' : [0, '', 'none://Section.py' ]
           }


import wx
import os
import sys
import string
from lxml import etree
from Section import *

def GetRelative(registered, segmented):
    if len(registered) != len(segmented):
        raise ValueError("both sets shoudl havethe same size")
    result = []
    
    for i in range(len(registered)):
        delta = registered[i] - segmented[i]
        result.append( ( 1.0 + delta * delta) / ( 1.0 + segmented[i] * segmented[i])

    return result



app = wx.App()
file = open(sys.argv[1], "r")
strset = file.read()
root = etree.XML(strset)
file.close()
set = WorkSet(node=root)

slice_nr = string.atoi(sys.argv[2])

if slice_nr >=  set.GetFrameNumber():
    print "Frame number  %d not from range:[0-%d)" % (slice_nr, set.GetFrameNumber())
    exit(1)

if len(sys.argv) > 3:
    zoom = atof(sys.argv[3])
else:
    zoom = 1.0

if len(sys.argv) > 4:
    shift = atof(sys.argv[4])
else:
    shift = 0.0


masks = []
slice = set.GetFrame(slice_nr)
image = slice.GetImageZoomed(zoom)
size = image.GetSize()


for sec in slice.Sections:
    masks.append(sec.GetMask(size, zoom))
    
    
for i in range(set.GetFrameNumber()):
    frame = set.GetFrame(i)
    image = frame.GetImageZoomed(zoom)
    p_ref = frame.GetMaskedIntensities(zoom)
    
    if len(p_ref) != 6:
        print "slice %d not segmented like it should be" % i
        exit(2) 
    
    p = []
    for mask in masks:
        val = 0.0
        n = len(mask)
        for m in mask:
            val = val + image.GetRed(m[0], m[1])
        p.append(val / n - shift)
    
    res = GetRelative(p, p_ref):

    print "%f %f %f %f %f %f" % (res[0], res[1], res[2], res[3], res[4], res[5])