Re: [Numpy-discussion] create a numpy array of images

2011-01-31 Thread totonixs...@gmail.com
I've been done that but with CT and MRI dicom files, and the cool
thing is that with numpy I can do something like this:

# getting axial slice
axial = slices[n,:,:]

# getting coronal slice
coronal = slices[:, n, :]

# getting sagital slice
sagital = slices[:,:, n]

On Sun, Jan 30, 2011 at 5:29 PM, Friedrich Romstedt
friedrichromst...@gmail.com wrote:
 2011/1/28 Christopher Barker chris.bar...@noaa.gov:
 On 1/28/11 7:01 AM, Asmi Shah wrote:
 I am using python for a while now and I have a requirement of creating a
 numpy array of microscopic tiff images ( this data is 3d, meaning there are
 100 z slices of 512 X 512 pixels.) How can I create an array of images?

 It's quite straightforward to create a 3-d array to hold this kind of data:

 image_block = np.empty((100, 512, 512), dtype=??)

 now you can load it up by using some lib (PIL, or ???) to load the tif
 images, and then:

 for i in images:
     image_block[i,:,:] = i

 Notice that since PIL 1.1.6, PIL Image objects support the numpy
 interface: http://effbot.org/zone/pil-changes-116.htm

 import PIL.Image
 im = PIL.Image.open('P1010102.JPG')
 im
 PIL.JpegImagePlugin.JpegImageFile image mode=RGB size=3264x2448 at 0x4CA0A8
 a = numpy.asarray(im)
 a.shape
 (2448, 3264, 3)
 a.dtype
 dtype('uint8')

 You can use the image just as any other ndarray:

 stack = numpy.empty((5, 2488, 3264, 3))
 stack[0] = im
 and so on

 for 5 images in a stack, notice that the dtype of the initially empty
 ndarray is float!

 It works also vice-versa:

 im_copy = PIL.Image.fromarray(a)

 but this seems to require integer-valued ndarrays as input, except
 when the ndarray is monochrome.

 This might be even simpler than the dtype proposed by Christopher.

 For more info on PIL: http://www.pythonware.com/library/pil/handbook/

 Friedrich
 ___
 NumPy-Discussion mailing list
 NumPy-Discussion@scipy.org
 http://mail.scipy.org/mailman/listinfo/numpy-discussion
___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Generator arrays

2011-01-31 Thread Neal Becker
I'm not sure how it applies to this discussion, but I'd just like to mention 
that a lot of interest (in c++ and d communities) has moved away from using 
iterators as the fundamental interface to containers and to ranges as the 
interface.

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


[Numpy-discussion] using loadtxt() for given number of rows?

2011-01-31 Thread Robert Cimrman
Hi,

I work with text files which contain several arrays separated by a few 
lines of other information, for example:

POINTS 4 float
-5.00e-01 -5.00e-01 0.00e+00
5.00e-01 -5.00e-01 0.00e+00
5.00e-01 5.00e-01 0.00e+00
-5.00e-01 5.00e-01 0.00e+00

CELLS 2 8
3 0 1 2
3 2 3 0

(yes, that's the legacy VTK format, but take it just as an example)

I have used custom Python code with loops to read similar files, so the 
speed was not too good. Now I wonder if it would be possible to use the 
numpy.loadtxt() function for the array-like parts. It supports passing 
an open file object in, which is good, but it wants to read the entire 
file, which does not work in this case.

It seems to me, that an additional parameter to loadtxt(), say nrows or 
numrows, would do the job, so that the function does not try reading the 
entire file. Another possibility would be to raise an exception as it 
is now, but also to return the data succesfully read so far.

What do you think? Is this worth a ticket?

r.
___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


[Numpy-discussion] Vectorize or rewrite function to work with array inputs?

2011-01-31 Thread DParker
I have several functions like the example below that I would like to make 
compatible with array inputs. The problem is the conditional statements 
give a ValueError: The truth value of an array with more than one element 
is ambiguous. Use a.any() or a.all(). I can use numpy.vectorize, but if 
possible I'd prefer to rewrite the function. Does anyone have any advice 
the best way to modify the code to accept array inputs? Thanks in advance 
for any assistance.


NAN = float('nan')

def air_gamma(t, far=0.0):

Specific heat ratio (gamma) of Air/JP8
t - static temperature, Rankine
[far] - fuel air ratio [- defaults to 0.0 (dry air)]
air_gamma - specific heat ratio

if far  0.:
return NAN
elif far  0.005:
if t  379. or t  4731.:
return NAN
else:
air_gamma = -3.472487e-22 * t ** 6. + 6.218811e-18 * t ** 5. - 
4.428098e-14 * t ** 4. + 1.569889e-10 * t ** 3. - 0.002753524 * t ** 
2. + 0.0001684666 * t + 1.368652
elif far  0.069:
if t  699. or t  4731.:
return NAN
else:
a6 = 4.114808e-20 * far ** 3. - 1.644588e-20 * far ** 2. + 
3.103507e-21 * far - 3.391308e-22
a5 = -6.819015e-16 * far ** 3. + 2.773945e-16 * far ** 2. - 
5.469399e-17 * far + 6.058125e-18
a4 = 4.684637e-12 * far ** 3. - 1.887227e-12 * far ** 2. + 
3.865306e-13 * far - 4.302534e-14
a3 = -0.0001700602 * far ** 3. + 0.6593809 * far 
** 2. - 0.1392629 * far + 1.520583e-10
a2 = 0.3431136 * far ** 3. - 0.1248285 * far ** 2. + 
0.02688007 * far - 0.002651616
a1 = -0.03792449 * far ** 3. + 0.01261025 * far ** 2. - 
0.002676877 * far + 0.0001580424
a0 = 13.65379 * far ** 3. - 3.311225 * far ** 2. + 0.3573201 * 
far + 1.372714
air_gamma = a6 * t ** 6. + a5 * t ** 5. + a4 * t ** 4. + a3 * 
t ** 3. + a2 * t ** 2. + a1 * t + a0
elif far = 0.069:
return NAN
else:
return NAN
return air_gamma

David Parker 
Chromalloy - TDAG___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Vectorize or rewrite function to work with array inputs?

2011-01-31 Thread Sebastian Berg
Hello,

On Mon, 2011-01-31 at 10:15 -0500, dpar...@chromalloy.com wrote:
 I have several functions like the example below that I would like to
 make compatible with array inputs. The problem is the conditional
 statements give a ValueError: The truth value of an array with more
 than one element is ambiguous. Use a.any() or a.all(). I can use
 numpy.vectorize, but if possible I'd prefer to rewrite the function.
 Does anyone have any advice the best way to modify the code to accept
 array inputs? Thanks in advance for any assistance. 
 

You can use binary indexing instead of the if:
condition = far  0
values[condition] = np.nan # set to NaN wherever far  0 is True.

or if you like I suppose you could put it into cython, add some typing
to avoid creating those binary arrays etc. all over to speed things up
more.

Regards,

Sebastian

 
 NAN = float('nan') 
 
 def air_gamma(t, far=0.0): 
  
 Specific heat ratio (gamma) of Air/JP8 
 t - static temperature, Rankine 
 [far] - fuel air ratio [- defaults to 0.0 (dry air)] 
 air_gamma - specific heat ratio 
  
 if far  0.: 
 return NAN 
 elif far  0.005: 
 if t  379. or t  4731.: 
 return NAN 
 else: 
 air_gamma = -3.472487e-22 * t ** 6. + 6.218811e-18 * t **
 5. - 4.428098e-14 * t ** 4. + 1.569889e-10 * t ** 3. - 0.002753524
 * t ** 2. + 0.0001684666 * t + 1.368652 
 elif far  0.069: 
 if t  699. or t  4731.: 
 return NAN 
 else: 
 a6 = 4.114808e-20 * far ** 3. - 1.644588e-20 * far ** 2. +
 3.103507e-21 * far - 3.391308e-22 
 a5 = -6.819015e-16 * far ** 3. + 2.773945e-16 * far ** 2.
 - 5.469399e-17 * far + 6.058125e-18 
 a4 = 4.684637e-12 * far ** 3. - 1.887227e-12 * far ** 2. +
 3.865306e-13 * far - 4.302534e-14 
 a3 = -0.0001700602 * far ** 3. + 0.6593809 *
 far ** 2. - 0.1392629 * far + 1.520583e-10 
 a2 = 0.3431136 * far ** 3. - 0.1248285 * far ** 2.
 + 0.02688007 * far - 0.002651616 
 a1 = -0.03792449 * far ** 3. + 0.01261025 * far ** 2. -
 0.002676877 * far + 0.0001580424 
 a0 = 13.65379 * far ** 3. - 3.311225 * far ** 2. +
 0.3573201 * far + 1.372714 
 air_gamma = a6 * t ** 6. + a5 * t ** 5. + a4 * t ** 4. +
 a3 * t ** 3. + a2 * t ** 2. + a1 * t + a0 
 elif far = 0.069: 
 return NAN 
 else: 
 return NAN 
 return air_gamma 
 
 David Parker 
 Chromalloy - TDAG
 ___
 NumPy-Discussion mailing list
 NumPy-Discussion@scipy.org
 http://mail.scipy.org/mailman/listinfo/numpy-discussion


___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Inversion of near singular matrices.

2011-01-31 Thread Sturla Molden
Den 31.01.2011 03:05, skrev Algis Kabaila:
 Actually, the structural engineer
 has no interest in trying to invert a singular matrix. However
 he/she is interested (or should be interested :)  )  when the
 square response matrix might approach singularity for this would
 signal instability.

I am sorry for having confused the issue by mentioning statistics. The 
mathematics (linear algebra) is of course the same. A singular matrix 
cannot be inverted by definition. The methods mentioned (SVD, Tikohonov 
regularization), as well as the transforms mentioned by Paul, will let 
you avoid numerical instability when matrices approach singularity 
(i.e. are very ill-conditioned).

OT: I think I know what structural engineering is. Back in 1994 I had to 
take a class in statikk (not sure what that translates to in English), 
with a textbook by Fritjof Irgens. From what I remember we did vector 
calculus to ensure the forces in a construction summed to 0, so that 
Newton's first law of motion would apply. It's unhealthy to be inside a 
building otherwise ;-)

Sturla Molden






___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Vectorize or rewrite function to work with array inputs?

2011-01-31 Thread eat
Hi,

On Mon, Jan 31, 2011 at 5:15 PM, dpar...@chromalloy.com wrote:

 I have several functions like the example below that I would like to make
 compatible with array inputs. The problem is the conditional statements give
 a *ValueError: The truth value of an array with more than one element is
 ambiguous. Use a.any() or a.all()*. I can use numpy.vectorize, but if
 possible I'd prefer to rewrite the function. Does anyone have any advice the
 best way to modify the code to accept array inputs? Thanks in advance for
 any assistance.


If I understod your question correctly, then air_gamma could be coded as:
def air_gamma_0(t, far=0.0):

Specific heat ratio (gamma) of Air/JP8
t - static temperature, Rankine
[far] - fuel air ratio [- defaults to 0.0 (dry air)]
air_gamma - specific heat ratio

if far 0.:
return NAN
elif far  0.005:
ag= air_gamma_1(t)
ag[np.logical_or(t 379., t 4731.)]= NAN
return ag
elif far 0.069:
ag= air_gamma_2(t, far)
ag[np.logical_or(t 699., t 4731.)]= NAN
return ag
else:
return NAN
Rest of the code is in the attachment.


My two cents,
eat




 NAN = float('nan')

 def air_gamma(t, far=0.0):
 
 Specific heat ratio (gamma) of Air/JP8
 t - static temperature, Rankine
 [far] - fuel air ratio [- defaults to 0.0 (dry air)]
 air_gamma - specific heat ratio
 
 if far  0.:
 return NAN
 elif far  0.005:
 if t  379. or t  4731.:
 return NAN
 else:
 air_gamma = -3.472487e-22 * t ** 6. + 6.218811e-18 * t ** 5. -
 4.428098e-14 * t ** 4. + 1.569889e-10 * t ** 3. - 0.002753524 * t ** 2.
 + 0.0001684666 * t + 1.368652
 elif far  0.069:
 if t  699. or t  4731.:
 return NAN
 else:
 a6 = 4.114808e-20 * far ** 3. - 1.644588e-20 * far ** 2. +
 3.103507e-21 * far - 3.391308e-22
 a5 = -6.819015e-16 * far ** 3. + 2.773945e-16 * far ** 2. -
 5.469399e-17 * far + 6.058125e-18
 a4 = 4.684637e-12 * far ** 3. - 1.887227e-12 * far ** 2. +
 3.865306e-13 * far - 4.302534e-14
 a3 = -0.0001700602 * far ** 3. + 0.6593809 * far **
 2. - 0.1392629 * far + 1.520583e-10
 a2 = 0.3431136 * far ** 3. - 0.1248285 * far ** 2. +
 0.02688007 * far - 0.002651616
 a1 = -0.03792449 * far ** 3. + 0.01261025 * far ** 2. -
 0.002676877 * far + 0.0001580424
 a0 = 13.65379 * far ** 3. - 3.311225 * far ** 2. + 0.3573201 *
 far + 1.372714
 air_gamma = a6 * t ** 6. + a5 * t ** 5. + a4 * t ** 4. + a3 * t
 ** 3. + a2 * t ** 2. + a1 * t + a0
 elif far = 0.069:
 return NAN
 else:
 return NAN
 return air_gamma

 David Parker
 Chromalloy - TDAG
 ___
 NumPy-Discussion mailing list
 NumPy-Discussion@scipy.org
 http://mail.scipy.org/mailman/listinfo/numpy-discussion




air_gamma.py
Description: Binary data
___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] numpy.linalg.svd documentation

2011-01-31 Thread Sturla Molden
Den 30.01.2011 21:40, skrev Charles R Harris:
 Well, strictly speaking, both documentations say the same thing, but 
 the old version was somewhat obfuscated. Either svd returns v.H and A 
 = dot(u*d, v.H) or svd returns v and A = dot(u*d,v). I think the 
 second is a clearer statement of the return value and the resulting 
 factorization, but I suppose some may hold a different opinion.

I agree that it is a clearer statement, since there is a difference 
between A = dot(u*d, v.H)  and A = dot(u*d, vH), and we actually have 
the latter. Still, the common definition of SVD is

 u * s * v.H = x

and not

 u * s * v = x

This might be a bit confusing for those expecting the conjugate 
transpose of v in the decomposition.

Be aware that Matlab's SVD is

[u,s,v] = svd(x)

so that u*s*v' = x.

Clearly Matlab and NumPy differ on the definition of v here, with Matlab 
following the common convention. That is why I prefer the old notation

u, s, vH = np.linalg.svd(x)

v = vH.H

as it leaves no room for confusion.

As long as the behaviour of SVD has not changed, none of my SVD code 
will break. That was what worried me most :-)

Sturla








___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Inversion of near singular matrices.

2011-01-31 Thread Algis Kabaila
On Tuesday 01 February 2011 03:27:22 Sturla Molden wrote:
 Den 31.01.2011 03:05, skrev Algis Kabaila:
  Actually, the structural engineer
  has no interest in trying to invert a singular matrix.
  However he/she is interested (or should be interested :) 
  )  when the square response matrix might approach
  singularity for this would signal instability.
 
 I am sorry for having confused the issue by mentioning
 statistics. The mathematics (linear algebra) is of course
 the same. A singular matrix cannot be inverted by
 definition. The methods mentioned (SVD, Tikohonov
 regularization), as well as the transforms mentioned by
 Paul, will let you avoid numerical instability when matrices
 approach singularity (i.e. are very ill-conditioned).
 
 OT: I think I know what structural engineering is. Back in
 1994 I had to take a class in statikk (not sure what that
 translates to in English), with a textbook by Fritjof
 Irgens. From what I remember we did vector calculus to
 ensure the forces in a construction summed to 0, so that
 Newton's first law of motion would apply. It's unhealthy to
 be inside a building otherwise ;-)
 
 Sturla Molden
 
I would guess that statikk is statics, the subject of 
conditions of equilibrium.  

Yes, teaching is not for the faint hearted... Particularly in 
foreign areas.  Just to put your mind at ese - it is important 
to have some idea of statistics even in simplest engineering 
structures, such as those made up of statically determinate 
trusses.  (A truss is made up of members that are pin jointed, 
or are imagined to be pin jointed.  Because of the pin joints, 
each member can only be subjected to an axial force. My next 
code snippet will show the vagaries of analisis of statically 
determinate trusses).

Before I can really ask my next question, I should know what 
matrix norms are used for the calculation of matrix condition 
number in numpy.linalg.  You see, I tried to compare it with a 
condition number found in an undergraduate text book and got a 
totally different number. 

So if you know that and are able to explain it in simple terms 
so that even engineers can understand it, it will be greatly 
appreciated. 

Al.

-- 
Algis
http://akabaila.pcug.org.au/StructuralAnalysis.pdf
___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] create a numpy array of images

2011-01-31 Thread Zachary Pincus
 I am using python for a while now and I have a requirement of  
 creating a
 numpy array of microscopic tiff images ( this data is 3d, meaning  
 there are
 100 z slices of 512 X 512 pixels.) How can I create an array of  
 images?

 It's quite straightforward to create a 3-d array to hold this kind  
 of data:

 image_block = np.empty((100, 512, 512), dtype=??)

 now you can load it up by using some lib (PIL, or ???) to load the  
 tif
 images, and then:

 for i in images:
 image_block[i,:,:] = i

 Notice that since PIL 1.1.6, PIL Image objects support the numpy
 interface: http://effbot.org/zone/pil-changes-116.htm

For even longer than this, PIL has been somewhat broken with regard to  
16-bit images (very common in microscopy); you may run into strange  
byte-ordering issues that scramble the data on reading or writing.  
Also, PIL's numpy interface is somewhat broken in similar ways.  
(Numerous people have provided patches to PIL, but these are never  
incorporated into any releases, as far as I can tell.)

So try PIL, but if the images come out all wrong, you might want to  
check out the scikits.image package, which has hooks for various other  
image read/write tools.

Zach
___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] using loadtxt() for given number of rows?

2011-01-31 Thread Christopher Barker
On 1/31/11 4:39 AM, Robert Cimrman wrote:
 I work with text files which contain several arrays separated by a few
 lines of other information, for example:

 POINTS 4 float
 -5.00e-01 -5.00e-01 0.00e+00
 5.00e-01 -5.00e-01 0.00e+00
 5.00e-01 5.00e-01 0.00e+00
 -5.00e-01 5.00e-01 0.00e+00

 CELLS 2 8
 3 0 1 2
 3 2 3 0

 I have used custom Python code with loops to read similar files, so the
 speed was not too good. Now I wonder if it would be possible to use the
 numpy.loadtxt() function for the array-like parts. It supports passing
 an open file object in, which is good, but it wants to read the entire
 file, which does not work in this case.

 It seems to me, that an additional parameter to loadtxt(), say nrows or
 numrows, would do the job,

I agree that that would be a useful feature. However, I'm not sure it 
would help performance much -- I think loadtxt is written in python as well.

One option in the meantime. If you know how many rows, you presumable 
know how many items on each row. IN that case, you can use:

np.fromfile(open_file, sep=' ', count=num_items_to_read)

It'll only work for multi-line text if the separator is whitespace, 
which it was in your example. But if it does, it should be pretty fast.

-Chris


-- 
Christopher Barker, Ph.D.
Oceanographer

Emergency Response Division
NOAA/NOS/ORR(206) 526-6959   voice
7600 Sand Point Way NE   (206) 526-6329   fax
Seattle, WA  98115   (206) 526-6317   main reception

chris.bar...@noaa.gov
___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Strange behaviour of numpy.asarray() in corner case

2011-01-31 Thread Friedrich Romstedt
2011/1/28 Friedrich Romstedt friedrichromst...@gmail.com:
 numpy.asarray([X(), numpy.asarray([1, 1])]).shape
 (2,)
 numpy.asarray([numpy.asarray([1, 1]), X()]).shape
 ()

Does noone have an opinion about this?  Shall I file a ticket?

Friedrich
___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Strange behaviour of numpy.asarray() in corner case

2011-01-31 Thread Warren Weckesser
On Mon, Jan 31, 2011 at 8:32 PM, Friedrich Romstedt 
friedrichromst...@gmail.com wrote:

 2011/1/28 Friedrich Romstedt friedrichromst...@gmail.com:
  numpy.asarray([X(), numpy.asarray([1, 1])]).shape
  (2,)
  numpy.asarray([numpy.asarray([1, 1]), X()]).shape
  ()

 Does noone have an opinion about this?



Looks wrong  to me.


 Shall I file a ticket?



Yes.


Warren




 Friedrich
 ___
 NumPy-Discussion mailing list
 NumPy-Discussion@scipy.org
 http://mail.scipy.org/mailman/listinfo/numpy-discussion

___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion


Re: [Numpy-discussion] Strange behaviour of numpy.asarray() in corner case

2011-01-31 Thread Friedrich Romstedt
2011/2/1 Warren Weckesser warren.weckes...@enthought.com:
  Shall I file a ticket?

 Yes.

Ok, #1730: http://projects.scipy.org/numpy/ticket/1730.

Thanks,
Friedrich
___
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion