Re: [Numpy-discussion] create a numpy array of images
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
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?
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?
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?
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.
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?
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
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.
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
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?
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/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
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/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