Search Results

Search found 236 results on 10 pages for 'p np'.

Page 2/10 | < Previous Page | 1 2 3 4 5 6 7 8 9 10  | Next Page >

  • Why do I get rows of zeros in my 2D fft?

    - by Nicholas Pringle
    I am trying to replicate the results from a paper. "Two-dimensional Fourier Transform (2D-FT) in space and time along sections of constant latitude (east-west) and longitude (north-south) were used to characterize the spectrum of the simulated flux variability south of 40degS." - Lenton et al(2006) The figures published show "the log of the variance of the 2D-FT". I have tried to create an array consisting of the seasonal cycle of similar data as well as the noise. I have defined the noise as the original array minus the signal array. Here is the code that I used to plot the 2D-FT of the signal array averaged in latitude: import numpy as np from numpy import ma from matplotlib import pyplot as plt from Scientific.IO.NetCDF import NetCDFFile ### input directory indir = '/home/nicholas/data/' ### get the flux data which is in ### [time(5day ave for 10 years),latitude,longitude] nc = NetCDFFile(indir + 'CFLX_2000_2009.nc','r') cflux_southern_ocean = nc.variables['Cflx'][:,10:50,:] cflux_southern_ocean = ma.masked_values(cflux_southern_ocean,1e+20) # mask land nc.close() cflux = cflux_southern_ocean*1e08 # change units of data from mmol/m^2/s ### create an array that consists of the seasonal signal fro each pixel year_stack = np.split(cflux, 10, axis=0) year_stack = np.array(year_stack) signal_array = np.tile(np.mean(year_stack, axis=0), (10, 1, 1)) signal_array = ma.masked_where(signal_array > 1e20, signal_array) # need to mask ### average the array over latitude(or longitude) signal_time_lon = ma.mean(signal_array, axis=1) ### do a 2D Fourier Transform of the time/space image ft = np.fft.fft2(signal_time_lon) mgft = np.abs(ft) ps = mgft**2 log_ps = np.log(mgft) log_mgft= np.log(mgft) Every second row of the ft consists completely of zeros. Why is this? Would it be acceptable to add a randomly small number to the signal to avoid this. signal_time_lon = signal_time_lon + np.random.randint(0,9,size=(730, 182))*1e-05 EDIT: Adding images and clarify meaning The output of rfft2 still appears to be a complex array. Using fftshift shifts the edges of the image to the centre; I still have a power spectrum regardless. I expect that the reason that I get rows of zeros is that I have re-created the timeseries for each pixel. The ft[0, 0] pixel contains the mean of the signal. So the ft[1, 0] corresponds to a sinusoid with one cycle over the entire signal in the rows of the starting image. Here are is the starting image using following code: plt.pcolormesh(signal_time_lon); plt.colorbar(); plt.axis('tight') Here is result using following code: ft = np.fft.rfft2(signal_time_lon) mgft = np.abs(ft) ps = mgft**2 log_ps = np.log1p(mgft) plt.pcolormesh(log_ps); plt.colorbar(); plt.axis('tight') It may not be clear in the image but it is only every second row that contains completely zeros. Every tenth pixel (log_ps[10, 0]) is a high value. The other pixels (log_ps[2, 0], log_ps[4, 0] etc) have very low values.

    Read the article

  • 3 dimensional bin packing algorithms

    - by BuschnicK
    I'm faced with a 3 dimensional bin packing problem and am currently conducting some preliminary research as to which algorithms/heuristics are currently yielding the best results. Since the problem is NP hard I do not expect to find the optimal solution in every case, but I was wondering: 1) what are the best exact solvers? Branch and Bound? What problem instance sizes can I expect to solve with reasonable computing resources? 2) what are the best heuristic solvers? 3) What off-the-shelf solutions exist to conduct some experiments with?

    Read the article

  • numpy calling sse2 via ctypes

    - by Daniel
    Hello, In brief, I am trying to call into a shared library from python, more specifically, from numpy. The shared library is implemented in C using sse2 instructions. Enabling optimisation, i.e. building the library with -O2 or –O1, I am facing strange segfaults when calling into the shared library via ctypes. Disabling optimisation (-O0), everything works out as expected, as is the case when linking the library to a c-program directly (optimised or not). Attached you find a snipped which exhibits the delineated behaviour on my system. With optimisation enabled, gdb reports a segfault in __builtin_ia32_loadupd (__P) at emmintrin.h:113. The value of __P is reported as optimised out. test.c: #include <emmintrin.h> #include <complex.h> void test(const int m, const double* x, double complex* y) { int i; __m128d _f, _x, _b; double complex f __attribute__( (aligned(16)) ); double complex b __attribute__( (aligned(16)) ); __m128d* _p; b = 1; _b = _mm_loadu_pd( (double *) &b ); _p = (__m128d*) y; for(i=0; i<m; ++i) { f = cexp(-I*x[i]); _f = _mm_loadu_pd( (double *) &f ); _x = _mm_loadu_pd( (double *) &x[i] ); _f = _mm_shuffle_pd(_f, _f, 1); *_p = _mm_add_pd(*_p, _f); *_p = _mm_add_pd(*_p, _x); *_p = _mm_mul_pd(*_p,_b); _p++; } return; } Compiler flags: gcc -o libtest.so -shared -std=c99 -msse2 -fPIC -O2 -g -lm test.c test.py: import numpy as np import os def zerovec_aligned(nr, dtype=np.float64, boundary=16): '''Create an aligned array of zeros. ''' size = nr * np.dtype(dtype).itemsize tmp = np.zeros(size + boundary, dtype=np.uint8) address = tmp.__array_interface__['data'][0] offset = boundary - address % boundary return tmp[offset:offset + size].view(dtype=dtype) lib = np.ctypeslib.load_library('libtest', '.' ) lib.test.restype = None lib.test.argtypes = [np.ctypeslib.ctypes.c_int, np.ctypeslib.ndpointer(np.float64, flags=('C', 'A') ), np.ctypeslib.ndpointer(np.complex128, flags=('C', 'A', 'W') )] n = 13 y = zerovec_aligned(n, dtype=np.complex128) x = np.ones(n, dtype=np.float64) # x = zerovec_aligned(n, dtype=np.float64) # x[:] = 1. lib.test(n,x,y) My system: Ubuntu Linux i686 2.6.31-22-generic Compiler: gcc (Ubuntu 4.4.1-4ubuntu9) Python: Python 2.6.4 (r264:75706, Dec 7 2009, 18:45:15) [GCC 4.4.1] Numpy: 1.4.0 I have taken provisions (cf. python code) that y is aligned and the alignment of x should not matter (I think; explicitly aligning x does not solve the problem though). Note also that i use _mm_loadu_pd instead of _mm_load_pd when loading b and f. For the C-only version _mm_load_pd works (as expected). However, when calling the function via ctypes using _mm_load_pd always segfaults (independent of optimisation). I have tried several days to sort out this issue without success ... and I am on the verge beating my monitor to death. Any input welcome. Daniel

    Read the article

  • Minimum cost strongly connected digraph

    - by Kazoom
    I have a digraph which is strongly connected (i.e. there is a path from i to j and j to i for each pair of nodes (i, j) in the graph G). I wish to find a strongly connected graph out of this graph such that the sum of all edges is the least. To put it differently, I need to get rid of edges in such a way that after removing them, the graph will still be strongly connected and of least cost for the sum of edges. I think it's an NP hard problem. I'm looking for an optimal solution, not approximation, for a small set of data like 20 nodes. Edit A more general description: Given a grap G(V,E) find a graph G'(V,E') such that if there exists a path from v1 to v2 in G than there also exists a path between v1 and v2 in G' and sum of each ei in E' is the least possible. so its similar to finding a minimum equivalent graph, only here we want to minimize the sum of edge weights rather than sum of edges. Edit: My approach so far: I thought of solving it using TSP with multiple visits, but it is not correct. My goal here is to cover each city but using a minimum cost path. So, it's more like the cover set problem, I guess, but I'm not exactly sure. I'm required to cover each and every city using paths whose total cost is minimum, so visiting already visited paths multiple times does not add to the cost.

    Read the article

  • Efficiently generate numpy array from list comprehension output?

    - by shootingstars
    Is there a more efficient way than using numpy.asarray() to generate an array from output in the form of a list? This appears to be copying everything in memory, which doesn't seem like it would be that efficient with very large arrays. (Updated) Example: import numpy as np a1 = np.array([1,2,3,4,5,6,7,8,9,10]) # pretend this has thousands of elements a2 = np.array([3,7,8]) results = np.asarray([np.amax(np.where(a1 > element)) for element in a2])

    Read the article

  • Traveling Salesman in polynomial time

    - by Andres
    Problem Description: Write a program in any language with the shortest number of characters that solves the following problem: Given a collection of cities and the cost of travel between each pair of them, find the cheapest way of visiting all of the cities and returning to your starting point. All solutions must take at most polynomial time. Input will be in the form of a text file named simply "i". The text file will contain data in the following format: city1# city2# cost$ cityA# cityB# cost$ ... Each element in the file is separated with a space, and there is a newline at the end of every line. Code count does not include whitespace. Solutions that take longer than polynomial time to find will not be accepted

    Read the article

  • Deterministic Annealing Code

    - by wade
    I would like to find an open source example of a code for deterministic annealing. It can be in almost any language: C, C++, MatLab/Octave, Fortran. I have already found a MatLab code for simulated annealing, so MatLab would be best. Here is a paper that describes the algorithm: http://www.google.com/url?sa=t&source=web&ct=res&cd=1&ved=0CB8QFjAA&url=http%3A%2F%2Fvandamteaching.googlepages.com%2FABriefIntroductionToDeterministicAnn.pdf&ei=DiLiS8qZFI7AMozB1JED&usg=AFQjCNHLps7HRWXLNN5rAX5aJ5BsJbcHuQ&sig2=YSokUTOs0UszAFZ9TDiJgQ

    Read the article

  • Optimizing a Parking Lot Problem. What algorithims should I use to fit the most amount of cars in th

    - by Adam Gent
    What algorithms (brute force or not) would I use to put in as many cars (assume all cars are the same size) in a parking lot so that there is at least one exit (from the container) and a car cannot be blocked. Or can someone show me an example of this problem solved programmatically. The parking lot varies in shape would be nice but if you want to assume its some invariant shape that is fine. Another Edit: Assume that driving distance in the parking lot is not a factor (although it would be totally awesome if it was weighted factor to number of cars in lot). Another Edit: Assume 2 Dimensional (no cranes or driving over cars). Another Edit: You cannot move cars around once they are parked (its not a valet parking lot). I hope the question is specific enough now.

    Read the article

  • Using scipy.interpolate.splrep function

    - by Koustav Ghosal
    I am trying to fit a cubic spline to a given set of points. My points are not ordered. I CANNOT sort or reorder the points, since I need that information. But since the function scipy.interpolate.splrep works only on non-duplicate and monotonically increasing points I have defined a function that maps the x-coordinates to a monotonically increasing space. My old points are: xpoints=[4913.0, 4912.0, 4914.0, 4913.0, 4913.0, 4913.0, 4914.0, 4915.0, 4918.0, 4921.0, 4925.0, 4932.0, 4938.0, 4945.0, 4950.0, 4954.0, 4955.0, 4957.0, 4956.0, 4953.0, 4949.0, 4943.0, 4933.0, 4921.0, 4911.0, 4898.0, 4886.0, 4874.0, 4865.0, 4858.0, 4853.0, 4849.0, 4848.0, 4849.0, 4851.0, 4858.0, 4864.0, 4869.0, 4877.0, 4884.0, 4893.0, 4903.0, 4913.0, 4923.0, 4935.0, 4947.0, 4959.0, 4970.0, 4981.0, 4991.0, 5000.0, 5005.0, 5010.0, 5015.0, 5019.0, 5020.0, 5021.0, 5023.0, 5025.0, 5027.0, 5027.0, 5028.0, 5028.0, 5030.0, 5031.0, 5033.0, 5035.0, 5037.0, 5040.0, 5043.0] ypoints=[10557.0, 10563.0, 10567.0, 10571.0, 10575.0, 10577.0, 10578.0, 10581.0, 10582.0, 10582.0, 10582.0, 10581.0, 10578.0, 10576.0, 10572.0, 10567.0, 10560.0, 10550.0, 10541.0, 10531.0, 10520.0, 10511.0, 10503.0, 10496.0, 10490.0, 10487.0, 10488.0, 10488.0, 10490.0, 10495.0, 10504.0, 10513.0, 10523.0, 10533.0, 10542.0, 10550.0, 10556.0, 10559.0, 10560.0, 10559.0, 10555.0, 10550.0, 10543.0, 10533.0, 10522.0, 10514.0, 10505.0, 10496.0, 10490.0, 10486.0, 10482.0, 10481.0, 10482.0, 10486.0, 10491.0, 10497.0, 10506.0, 10516.0, 10524.0, 10534.0, 10544.0, 10552.0, 10558.0, 10564.0, 10569.0, 10573.0, 10576.0, 10578.0, 10581.0, 10582.0] Plots: The code for the mapping function and interpolation is: xnew=[] ynew=ypoints for c3,i in enumerate(xpoints): if np.isfinite(np.log(i*pow(2,c3))): xnew.append(np.log(i*pow(2,c3))) else: if c==0: xnew.append(np.random.random_sample()) else: xnew.append(xnew[c3-1]+np.random.random_sample()) xnew=np.asarray(xnew) ynew=np.asarray(ynew) constant1=10.0 nknots=len(xnew)/constant1 idx_knots = (np.arange(1,len(xnew)-1,(len(xnew)-2)/np.double(nknots))).astype('int') knots = [xnew[i] for i in idx_knots] knots = np.asarray(knots) int_range=np.linspace(min(xnew),max(xnew),len(xnew)) tck = interpolate.splrep(xnew,ynew,k=3,task=-1,t=knots) y1= interpolate.splev(int_range,tck,der=0) The code is throwing an error at the function interpolate.splrep() for some set of points like the above one. The error is: File "/home/neeraj/Desktop/koustav/res/BOS5/fit_spline3.py", line 58, in save_spline_f tck = interpolate.splrep(xnew,ynew,k=3,task=-1,t=knots) File "/usr/lib/python2.7/dist-packages/scipy/interpolate/fitpack.py", line 465, in splrep raise _iermessier(_iermess[ier][0]) ValueError: Error on input data But for other set of points it works fine. For example for the following set of points. xpoints=[1629.0, 1629.0, 1629.0, 1629.0, 1629.0, 1629.0, 1629.0, 1629.0, 1629.0, 1629.0, 1629.0, 1629.0, 1629.0, 1629.0, 1629.0, 1629.0, 1630.0, 1630.0, 1630.0, 1631.0, 1631.0, 1631.0, 1631.0, 1630.0, 1629.0, 1629.0, 1629.0, 1628.0, 1627.0, 1627.0, 1625.0, 1624.0, 1624.0, 1623.0, 1620.0, 1618.0, 1617.0, 1616.0, 1615.0, 1614.0, 1614.0, 1612.0, 1612.0, 1612.0, 1611.0, 1610.0, 1609.0, 1608.0, 1607.0, 1607.0, 1603.0, 1602.0, 1602.0, 1601.0, 1601.0, 1600.0, 1599.0, 1598.0] ypoints=[10570.0, 10572.0, 10572.0, 10573.0, 10572.0, 10572.0, 10571.0, 10570.0, 10569.0, 10565.0, 10564.0, 10563.0, 10562.0, 10560.0, 10558.0, 10556.0, 10554.0, 10551.0, 10548.0, 10547.0, 10544.0, 10542.0, 10541.0, 10538.0, 10534.0, 10532.0, 10531.0, 10528.0, 10525.0, 10522.0, 10519.0, 10517.0, 10516.0, 10512.0, 10509.0, 10509.0, 10507.0, 10504.0, 10502.0, 10500.0, 10501.0, 10499.0, 10498.0, 10496.0, 10491.0, 10492.0, 10488.0, 10488.0, 10488.0, 10486.0, 10486.0, 10485.0, 10485.0, 10486.0, 10483.0, 10483.0, 10482.0, 10480.0] Plots: Can anybody suggest what's happening ?? Thanks in advance......

    Read the article

  • [numpy] storing record arrays in object arrays

    - by Peter Prettenhofer
    I'd like to convert a list of record arrays -- dtype is (uint32, float32) -- into a numpy array of dtype np.object: X = np.array(instances, dtype = np.object) where instances is a list of arrays with data type np.dtype([('f0', '<u4'), ('f1', '<f4')]). However, the above statement results in an array whose elements are also of type np.object: X[0] array([(67111L, 1.0), (104242L, 1.0)], dtype=object) Does anybody know why? The following statement should be equivalent to the above but gives the desired result: X = np.empty((len(instances),), dtype = np.object) X[:] = instances X[0] array([(67111L, 1.0), (104242L, 1.0), dtype=[('f0', '<u4'), ('f1', '<f4')]) thanks & best regards, peter

    Read the article

  • Matplotlib plotting non uniform data in 3D surface

    - by Raj Tendulkar
    I have a simple code to plot the points in 3D for Matplotlib as below - from mpl_toolkits.mplot3d import axes3d import matplotlib.pyplot as plt import numpy as np from numpy import genfromtxt import csv fig = plt.figure() ax = fig.add_subplot(111, projection='3d') my_data = genfromtxt('points1.csv', delimiter=',') points1X = my_data[:,0] points1Y = my_data[:,1] points1Z = my_data[:,2] ## I remove the header of the CSV File. points1X = np.delete(points1X, 0) points1Y = np.delete(points1Y, 0) points1Z = np.delete(points1Z, 0) # Convert the array to 1D array points1X = np.reshape(points1X,points1X.size) points1Y = np.reshape(points1Y,points1Y.size) points1Z = np.reshape(points1Z,points1Z.size) my_data = genfromtxt('points2.csv', delimiter=',') points2X = my_data[:,0] points2Y = my_data[:,1] points2Z = my_data[:,2] ## I remove the header of the CSV File. points2X = np.delete(points2X, 0) points2Y = np.delete(points2Y, 0) points2Z = np.delete(points2Z, 0) # Convert the array to 1D array points2X = np.reshape(points2X,points2X.size) points2Y = np.reshape(points2Y,points2Y.size) points2Z = np.reshape(points2Z,points2Z.size) ax.plot(points1X, points1Y, points1Z, 'd', markersize=8, markerfacecolor='red', label='points1') ax.plot(points2X, points2Y, points2Z, 'd', markersize=8, markerfacecolor='blue', label='points2') plt.show() My problem is that I tried to make a decent surface plot out of these data points that I have. I already tried to use ax.plot_surface() function to make it look nice. For this I eliminated some points, and recalculated the matrix kind of input needed by this function. However, the graph I generated was far more difficult to interpret and understand. So there might be 2 possibilities: either I am not using the function correctly, or otherwise, the data I am trying to plot is not good for the surface plot. What I was expecting was 3D graph which would have an effect similar to that we have of 3D pie chart. We see that one piece (that which is extracted out) is part of another piece. I was not expecting it to be exactly same like that, but some kind of effect like that. What I would like to ask is: Do you think it will be possible to make such 3D graph? Is there any way better, I could express my data in 3 dimension? Here are the 2 files - points1.csv Dim1,Dim2,Dim3 3,8,1 3,8,2 3,8,3 3,8,4 3,8,5 3,9,1 3,9,2 3,9,3 3,9,4 3,9,5 3,10,1 3,10,2 3,10,3 3,10,4 3,10,5 3,11,1 3,11,2 3,11,3 3,11,4 3,11,5 3,12,1 3,12,2 3,13,1 3,13,2 3,14,1 3,14,2 3,15,1 3,15,2 3,16,1 3,16,2 3,17,1 3,17,2 3,18,1 3,18,2 4,8,1 4,8,2 4,8,3 4,8,4 4,8,5 4,9,1 4,9,2 4,9,3 4,9,4 4,9,5 4,10,1 4,10,2 4,10,3 4,10,4 4,10,5 4,11,1 4,11,2 4,11,3 4,11,4 4,11,5 4,12,1 4,13,1 4,14,1 4,15,1 4,16,1 4,17,1 4,18,1 5,8,1 5,8,2 5,8,3 5,8,4 5,8,5 5,9,1 5,9,2 5,9,3 5,9,4 5,9,5 5,10,1 5,10,2 5,10,3 5,10,4 5,10,5 5,11,1 5,11,2 5,11,3 5,11,4 5,11,5 5,12,1 5,13,1 5,14,1 5,15,1 5,16,1 5,17,1 5,18,1 6,8,1 6,8,2 6,8,3 6,8,4 6,8,5 6,9,1 6,9,2 6,9,3 6,9,4 6,9,5 6,10,1 6,11,1 6,12,1 6,13,1 6,14,1 6,15,1 6,16,1 6,17,1 6,18,1 7,8,1 7,8,2 7,8,3 7,8,4 7,8,5 7,9,1 7,9,2 7,9,3 7,9,4 7,9,5 and points2.csv Dim1,Dim2,Dim3 3,12,3 3,12,4 3,12,5 3,13,3 3,13,4 3,13,5 3,14,3 3,14,4 3,14,5 3,15,3 3,15,4 3,15,5 3,16,3 3,16,4 3,16,5 3,17,3 3,17,4 3,17,5 3,18,3 3,18,4 3,18,5 4,12,2 4,12,3 4,12,4 4,12,5 4,13,2 4,13,3 4,13,4 4,13,5 4,14,2 4,14,3 4,14,4 4,14,5 4,15,2 4,15,3 4,15,4 4,15,5 4,16,2 4,16,3 4,16,4 4,16,5 4,17,2 4,17,3 4,17,4 4,17,5 4,18,2 4,18,3 4,18,4 4,18,5 5,12,2 5,12,3 5,12,4 5,12,5 5,13,2 5,13,3 5,13,4 5,13,5 5,14,2 5,14,3 5,14,4 5,14,5 5,15,2 5,15,3 5,15,4 5,15,5 5,16,2 5,16,3 5,16,4 5,16,5 5,17,2 5,17,3 5,17,4 5,17,5 5,18,2 5,18,3 5,18,4 5,18,5 6,10,2 6,10,3 6,10,4 6,10,5 6,11,2 6,11,3 6,11,4 6,11,5 6,12,2 6,12,3 6,12,4 6,12,5 6,13,2 6,13,3 6,13,4 6,13,5 6,14,2 6,14,3 6,14,4 6,14,5 6,15,2 6,15,3 6,15,4 6,15,5 6,16,2 6,16,3 6,16,4 6,16,5 6,17,2 6,17,3 6,17,4 6,17,5 6,18,2 6,18,3 6,18,4 6,18,5 7,10,1 7,10,2 7,10,3 7,10,4 7,10,5 7,11,1 7,11,2 7,11,3 7,11,4 7,11,5 7,12,1 7,12,2 7,12,3 7,12,4 7,12,5 7,13,1 7,13,2 7,13,3 7,13,4 7,13,5 7,14,1 7,14,2 7,14,3 7,14,4 7,14,5 7,15,1 7,15,2 7,15,3 7,15,4 7,15,5 7,16,1 7,16,2 7,16,3 7,16,4 7,16,5 7,17,1 7,17,2 7,17,3 7,17,4 7,17,5 7,18,1 7,18,2 7,18,3 7,18,4 7,18,5

    Read the article

  • Transferring data from 2d Dynamic array in C to CUDA and back

    - by Soumya
    I have a dynamically declared 2D array in my C program, the contents of which I want to transfer to a CUDA kernel for further processing. Once processed, I want to populate the dynamically declared 2D array in my C code with the CUDA processed data. I am able to do this with static 2D C arrays but not with dynamically declared C arrays. Any inputs would be welcome! I mean the dynamic array of dynamic arrays. The test code that I have written is as below. #include "cuda_runtime.h" #include "device_launch_parameters.h" #include <stdio.h> #include <conio.h> #include <math.h> #include <stdlib.h> const int nItt = 10; const int nP = 5; __device__ int d_nItt = 10; __device__ int d_nP = 5; __global__ void arr_chk(float *d_x_k, float *d_w_k, int row_num) { int index = (blockIdx.x * blockDim.x) + threadIdx.x; int index1 = (row_num * d_nP) + index; if ( (index1 >= row_num * d_nP) && (index1 < ((row_num +1)*d_nP))) //Modifying only one row data pertaining to one particular iteration { d_x_k[index1] = row_num * d_nP; d_w_k[index1] = index; } } float **mat_create2(int r, int c) { float **dynamicArray; dynamicArray = (float **) malloc (sizeof (float)*r); for(int i=0; i<r; i++) { dynamicArray[i] = (float *) malloc (sizeof (float)*c); for(int j= 0; j<c;j++) { dynamicArray[i][j] = 0; } } return dynamicArray; } /* Freeing memory - here only number of rows are passed*/ void cleanup2d(float **mat_arr, int x) { int i; for(i=0; i<x; i++) { free(mat_arr[i]); } free(mat_arr); } int main() { //float w_k[nItt][nP]; //Static array declaration - works! //float x_k[nItt][nP]; // if I uncomment this dynamic declaration and comment the static one, it does not work..... float **w_k = mat_create2(nItt,nP); float **x_k = mat_create2(nItt,nP); float *d_w_k, *d_x_k; // Device variables for w_k and x_k int nblocks, blocksize, nthreads; for(int i=0;i<nItt;i++) { for(int j=0;j<nP;j++) { x_k[i][j] = (nP*i); w_k[i][j] = j; } } for(int i=0;i<nItt;i++) { for(int j=0;j<nP;j++) { printf("x_k[%d][%d] = %f\t",i,j,x_k[i][j]); printf("w_k[%d][%d] = %f\n",i,j,w_k[i][j]); } } int size1 = nItt * nP * sizeof(float); printf("\nThe array size in memory bytes is: %d\n",size1); cudaMalloc( (void**)&d_x_k, size1 ); cudaMalloc( (void**)&d_w_k, size1 ); if((nP*nItt)<32) { blocksize = nP*nItt; nblocks = 1; } else { blocksize = 32; // Defines the number of threads running per block. Taken equal to warp size nthreads = blocksize; nblocks = ceil(float(nP*nItt) / nthreads); // Calculated total number of blocks thus required } for(int i = 0; i< nItt; i++) { cudaMemcpy( d_x_k, x_k, size1,cudaMemcpyHostToDevice ); //copy of x_k to device cudaMemcpy( d_w_k, w_k, size1,cudaMemcpyHostToDevice ); //copy of w_k to device arr_chk<<<nblocks, blocksize>>>(d_x_k,d_w_k,i); cudaMemcpy( x_k, d_x_k, size1, cudaMemcpyDeviceToHost ); cudaMemcpy( w_k, d_w_k, size1, cudaMemcpyDeviceToHost ); } printf("\nVerification after return from gpu\n"); for(int i = 0; i<nItt; i++) { for(int j=0;j<nP;j++) { printf("x_k[%d][%d] = %f\t",i,j,x_k[i][j]); printf("w_k[%d][%d] = %f\n",i,j,w_k[i][j]); } } cudaFree( d_x_k ); cudaFree( d_w_k ); cleanup2d(x_k,nItt); cleanup2d(w_k,nItt); getch(); return 0;

    Read the article

  • matplotlib and python multithread file processing

    - by Napseis
    I have a large number of files to process. I have written a script that get, sort and plot the datas I want. So far, so good. I have tested it and it gives the desired result. Then I wanted to do this using multithreading. I have looked into the doc and examples on the internet, and using one thread in my program works fine. But when I use more, at some point I get random matplotlib error, and I suspect some conflict there, even though I use a function with names for the plots, and iI can't see where the problem could be. Here is the whole script should you need more comment, i'll add them. Thank you. #!/usr/bin/python import matplotlib matplotlib.use('GTKAgg') import numpy as np from scipy.interpolate import griddata import matplotlib.pyplot as plt import matplotlib.colors as mcl from matplotlib import rc #for latex import time as tm import sys import threading import Queue #queue in 3.2 and Queue in 2.7 ! import pdb #the debugger rc('text', usetex=True)#for latex map=0 #initialize the map index. It will be use to index the array like this: array[map,[x,y]] time=np.zeros(1) #an array to store the time middle_h=np.zeros((0,3)) #x phi c #for the middle of the box current_file=open("single_void_cyl_periodic_phi_c_middle_h_out",'r') for line in current_file: if line.startswith('# === time'): map+=1 np.append(time,[float(line.strip('# === time '))]) elif line.startswith('#'): pass else: v=np.fromstring(line,dtype=float,sep=' ') middle_h=np.vstack( (middle_h,v[[1,3,4]]) ) current_file.close() middle_h=middle_h.reshape((map,-1,3)) #3d array: map, x, phi,c ##### def load_and_plot(): #will load a map file, and plot it along with the corresponding profile loaded before while not exit_flag: print("fecthing work ...") #try: if not tasks_queue.empty(): map_index=tasks_queue.get() print("----> working on map: %s" %map_index) x,y,zp=np.loadtxt("single_void_cyl_growth_periodic_post_map_"+str(map_index),unpack=True, usecols=[1, 2,3]) for i,el in enumerate(zp): if el<0.: zp[i]=0. xv=np.unique(x) yv=np.unique(y) X,Y= np.meshgrid(xv,yv) Z = griddata((x, y), zp, (X, Y),method='nearest') figure=plt.figure(num=map_index,figsize=(14, 8)) ax1=plt.subplot2grid((2,2),(0,0)) ax1.plot(middle_h[map_index,:,0],middle_h[map_index,:,1],'*b') ax1.grid(True) ax1.axis([-15, 15, 0, 1]) ax1.set_title('Profiles') ax1.set_ylabel(r'$\phi$') ax1.set_xlabel('x') ax2=plt.subplot2grid((2,2),(1,0)) ax2.plot(middle_h[map_index,:,0],middle_h[map_index,:,2],'*r') ax2.grid(True) ax2.axis([-15, 15, 0, 1]) ax2.set_ylabel('c') ax2.set_xlabel('x') ax3=plt.subplot2grid((2,2),(0,1),rowspan=2,aspect='equal') sub_contour=ax3.contourf(X,Y,Z,np.linspace(0,1,11),vmin=0.) figure.colorbar(sub_contour,ax=ax3) figure.savefig('single_void_cyl_'+str(map_index)+'.png') plt.close(map_index) tasks_queue.task_done() else: print("nothing left to do, other threads finishing,sleeping 2 seconds...") tm.sleep(2) # except: # print("failed this time: %s" %map_index+". Sleeping 2 seconds") # tm.sleep(2) ##### exit_flag=0 nb_threads=2 tasks_queue=Queue.Queue() threads_list=[] jobs=list(range(map)) #each job is composed of a map print("inserting jobs in the queue...") for job in jobs: tasks_queue.put(job) print("done") #launch the threads for i in range(nb_threads): working_bee=threading.Thread(target=load_and_plot) working_bee.daemon=True print("starting thread "+str(i)+' ...') threads_list.append(working_bee) working_bee.start() #wait for all tasks to be treated tasks_queue.join() #flip the flag, so the threads know it's time to stop exit_flag=1 for t in threads_list: print("waiting for threads %s to stop..."%t) t.join() print("all threads stopped")

    Read the article

  • How can I draw a log-normalized imshow plot with a colorbar representing the raw data in matplotlib

    - by Adam Fraser
    I'm using matplotlib to plot log-normalized images but I would like the original raw image data to be represented in the colorbar rather than the [0-1] interval. I get the feeling there's a more matplotlib'y way of doing this by using some sort of normalization object and not transforming the data beforehand... in any case, there could be negative values in the raw image. import matplotlib.pyplot as plt import numpy as np def log_transform(im): '''returns log(image) scaled to the interval [0,1]''' try: (min, max) = (im[im > 0].min(), im.max()) if (max > min) and (max > 0): return (np.log(im.clip(min, max)) - np.log(min)) / (np.log(max) - np.log(min)) except: pass return im a = np.ones((100,100)) for i in range(100): a[i] = i f = plt.figure() ax = f.add_subplot(111) res = ax.imshow(log_transform(a)) # the colorbar drawn shows [0-1], but I want to see [0-99] cb = f.colorbar(res) I've tried using cb.set_array, but that didn't appear to do anything, and cb.set_clim, but that rescales the colors completely. Thanks in advance for any help :)

    Read the article

  • How do I open a file with a program via a shortcut from the cmd prompt

    - by PassByReference0
    Here's my predicament: When I add a program's location to my PATH, I can do the following in cmd prompt to open a file in my current directory: notepad++ open_me.txt And this opens open_me.txt in notepad++. However, I don't want to have to add every single program I want to run to my path. What I want is to add a folder called C:\Users\Me\Documents\Programs to my path and just drop shortcuts to various programs into that folder and have them function the same as adding them to my path. So I dropped a link to notepad++.exe named "np" in my folder, and what I got was this: I have to run it with start np (instead of just np) But more importantly, if I try start np open_me.txt, it opens notepad++.exe but looks for open_me.txt in notepad++'s directory. How can I do this properly? (Also, I'd like to be opening notepad++.exe with the shortened name of np)

    Read the article

  • 2d convolution using python and numpy

    - by mikip
    Hi I am trying to perform a 2d convolution in python using numpy I have a 2d array as follows with kernel H_r for the rows and H_c for the columns data = np.zeros((nr, nc), dtype=np.float32) #fill array with some data here then convolve for r in range(nr): data[r,:] = np.convolve(data[r,:], H_r, 'same') for c in range(nc): data[:,c] = np.convolve(data[:,c], H_c, 'same') It does not produce the output that I was expecting, does this code look OK Thanks

    Read the article

  • Numpy/Python performing terribly vs. Matlab

    - by Nissl
    Novice programmer here. I'm writing a program that analyzes the relative spatial locations of points (cells). The program gets boundaries and cell type off an array with the x coordinate in column 1, y coordinate in column 2, and cell type in column 3. It then checks each cell for cell type and appropriate distance from the bounds. If it passes, it then calculates its distance from each other cell in the array and if the distance is within a specified analysis range it adds it to an output array at that distance. My cell marking program is in wxpython so I was hoping to develop this program in python as well and eventually stick it into the GUI. Unfortunately right now python takes ~20 seconds to run the core loop on my machine while MATLAB can do ~15 loops/second. Since I'm planning on doing 1000 loops (with a randomized comparison condition) on ~30 cases times several exploratory analysis types this is not a trivial difference. I tried running a profiler and array calls are 1/4 of the time, almost all of the rest is unspecified loop time. Here is the python code for the main loop: for basecell in range (0, cellnumber-1): if firstcelltype == np.array((cellrecord[basecell,2])): xloc=np.array((cellrecord[basecell,0])) yloc=np.array((cellrecord[basecell,1])) xedgedist=(xbound-xloc) yedgedist=(ybound-yloc) if xloc>excludedist and xedgedist>excludedist and yloc>excludedist and yedgedist>excludedist: for comparecell in range (0, cellnumber-1): if secondcelltype==np.array((cellrecord[comparecell,2])): xcomploc=np.array((cellrecord[comparecell,0])) ycomploc=np.array((cellrecord[comparecell,1])) dist=math.sqrt((xcomploc-xloc)**2+(ycomploc-yloc)**2) dist=round(dist) if dist>=1 and dist<=analysisdist: arraytarget=round(dist*analysisdist/intervalnumber) addone=np.array((spatialraw[arraytarget-1])) addone=addone+1 targetcell=arraytarget-1 np.put(spatialraw,[targetcell,targetcell],addone) Here is the matlab code for the main loop: for basecell = 1:cellnumber; if firstcelltype==cellrecord(basecell,3); xloc=cellrecord(basecell,1); yloc=cellrecord(basecell,2); xedgedist=(xbound-xloc); yedgedist=(ybound-yloc); if (xloc>excludedist) && (yloc>excludedist) && (xedgedist>excludedist) && (yedgedist>excludedist); for comparecell = 1:cellnumber; if secondcelltype==cellrecord(comparecell,3); xcomploc=cellrecord(comparecell,1); ycomploc=cellrecord(comparecell,2); dist=sqrt((xcomploc-xloc)^2+(ycomploc-yloc)^2); if (dist>=1) && (dist<=100.4999); arraytarget=round(dist*analysisdist/intervalnumber); spatialsum(1,arraytarget)=spatialsum(1,arraytarget)+1; end end end end end end Thanks!

    Read the article

  • Operand should contain 1 column(s) about insert into & select

    - by user1038890
    "insert into NodeProfileSections (profile_no, tpl_section_no) select (np.profile_no, tps.tpl_section_no) from NodeProfile np, TemplateProfileSection tps, TemplateProfile tp where np.hostname = '%s' AND np.role = '%s' AND tp.tpl_profile_no = tps.tpl_profile_no AND tp.tpl_name = '%s' AND tp.role = '%s' AND tps.tpl_section_name = '%s';" %(hostname, role, template_name, role, section_name) error_message = 'Operand should contain 1 column(s)' How to solve this problem?

    Read the article

  • Flex3 / Air 2: NativeProcess doesn't accepts standard input data (Error #2044 & #3218)

    - by Edward
    Hi: I'm trying to open cmd.exe on a new process and pass some code to programatically eject a device; but when trying to do this all I get is: "Error #2044: Unhandled IOErrorEvent:. text=Error #3218: Error while writing data to NativeProcess.standardInput." Here's my code: private var NP:NativeProcess = new NativeProcess(); private function EjectDevice():void { var RunDLL:File = new File("C:\\Windows\\System32\\cmd.exe"); var NPI:NativeProcessStartupInfo = new NativeProcessStartupInfo(); NPI.executable = RunDLL; NP.start(NPI); NP.addEventListener(Event.STANDARD_OUTPUT_CLOSE, CatchOutput, false, 0, true); NP.standardInput.writeUTFBytes("start C:\\Windows\\System32\\rundll32.exe shell32.dll,Control_RunDLL hotplug.dll"); NP.closeInput(); } I also tried with writeUTF instead of writeUTFBytes, but I still get the error. Does anyone have an idea of what I'm doing wrong?. Thanks for your time :) Edward.

    Read the article

  • Python Least-Squares Natural Splines

    - by Eldila
    I am trying to find a numerical package which will fit a natural which minimizes weighted least squares. There is a package in scipy which does what I want for unnatural splines. import numpy as np import matplotlib.pyplot as plt from scipy import interpolate import random x = np.arange(0,5,1.0/2) xs = np.arange(0,5,1.0/500) y = np.sin(x+1) for i in range(len(y)): y[i] += .2*random.random() - .1 knots = np.array([1,2,3,4]) tck = interpolate.splrep(x,y,s=1,k=3,t=knots,task=-1) ynew = interpolate.splev(xs,tck,der=0) plt.figure() plt.plot(xs,ynew,x,y,'x')

    Read the article

  • How to get parent node in Stanford's JavaNLP?

    - by roddik
    Hello. Suppose I have such chunk of a sentence: (NP (NP (DT A) (JJ single) (NN page)) (PP (IN in) (NP (DT a) (NN wiki) (NN website)))) At a certain moment of time I have a reference to (JJ single) and I want to get the NP node binding A single page. If I get it right, that NP is the parent of the node, A and page are its siblings and it has no children (?). When I try to use the .parent() method of a tree, I always get null. The API says that's because the implementation doesn't know how to determine the parent node. Another method of interest is .ancestor(int height, Tree root), but I don't know how to get the root of the node. In both cases, since the parser knows how to indent and group trees, it must know the "parent" tree, right? How can I get it? Thanks

    Read the article

  • subplot matplotlib wrong syntax

    - by madptr
    I am using matplotlib to subplot in a loop. For instance, i would like to subplot 49 data sets, and from the doc, i implemented it this way; import numpy as np import matplotlib.pyplot as plt X1=list(range(0,10000,1)) X1 = [ x/float(10) for x in X1 ] nb_mix = 2 parameters = [] for i in range(49): param = [] Y = [0] * len(X1) for j in range(nb_mix): mean = 5* (1 + (np.random.rand() * 2 - 1 ) * 0.5 ) var = 10* (1 + np.random.rand() * 2 - 1 ) scale = 5* ( 1 + (np.random.rand() * 2 - 1) * 0.5 ) Y = [ Y[k] + scale * np.exp(-((X1[k] - mean)/float(var))**2) for k in range(len(X1)) ] param = param + [[mean, var, scale]] ax = plt.subplot(7, 7, i + 1) ax.plot(X1, Y) parameters = parameters + [param] ax.show() However, i have an index out of range error from i=0 onwards. Where can i do better to have it works ? Thanks

    Read the article

< Previous Page | 1 2 3 4 5 6 7 8 9 10  | Next Page >