![(please configure the [header_logo] section in trac.ini)](/twice/chrome/site/your_project_logo.png){kind=logo}
- Timestamp:
- 08/22/14 15:06:02 (10 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/src/scripts_Laura/ARCTIC/Travail_CEN/cartesian_grid_test.py
r46 r56 32 32 # definition of the new cartesian grid (x, y) # 33 33 ############################################### 34 x0 = -3000. # min limit of grid35 x1 = 2500. # max limit of grid34 x0 = -3000. # min x limit of grid (in km) 35 x1 = 2500. # max x limit of grid (in km) 36 36 xvec = np.arange(x0, x1+dx, dx) 37 37 nx = len(xvec) 38 y0 = -3000. # min limit of grid39 y1 = 3000. # max limit of grid38 y0 = -3000. # min y limit of grid (in km) 39 y1 = 3000. # max y limit of grid (in km) 40 40 yvec = np.arange(y0, y1+dy, dy) 41 41 ny = len(yvec) 42 xgrid_cart, ygrid_cart= np.meshgrid(xvec, yvec) # new cartesian grid (centered on North pole )42 xgrid_cart, ygrid_cart= np.meshgrid(xvec, yvec) # new cartesian grid (centered on North pole [x, y] = [0, 0]) 43 43 44 44 … … 46 46 # calculation of the new gridding # 47 47 ################################### 48 zgrid_output = np.zeros([ny, nx, nb_days], float) 49 ngrid_output = np.zeros([ny, nx, nb_days], float) 50 z2grid_output = np.zeros([ny, nx, nb_days], float) 51 sigmagrid_output = np.zeros([ny, nx, nb_days], float) 52 bblat = nonzero(lati >= 65.) # only select high latitude values of z48 zgrid_output = np.zeros([ny, nx, nb_days], float) # 3D-array of data in new cartesian grid 49 ngrid_output = np.zeros([ny, nx, nb_days], float) # 3D-array of number of data in each grid cell of new cartesian grid 50 z2grid_output = np.zeros([ny, nx, nb_days], float) # 3D-array of data*data in new cartesian grid 51 sigmagrid_output = np.zeros([ny, nx, nb_days], float) # 3D-array of std of data in each grid cell of new cartesian grid 52 bblat = nonzero(lati >= 65.) # only select high latitude observations (above 65° lat) 53 53 for ijr in range (0, nb_days): # loop on time - here time = days 54 54 bbjour = nonzero(jour[bblat] == ijr + 1.)[0] … … 61 61 # origin of the (x, y) grid : (x=0, y=0) <=> (lon=0, lat=0) 62 62 L = len(z) 63 Rt = 6371. 64 alpha = (pi*Rt)/180. 63 Rt = 6371. # radius of earth 64 alpha = (pi*Rt)/180. # convertion from deg to rad angles 65 65 beta = pi/180. 66 x = np.zeros([L], float) 67 y = np.zeros([L], float) 66 x = np.zeros([L], float) # x coordinates of new cartesian grid 67 y = np.zeros([L], float) # y coordinates of new cartesian grid 68 68 for k in range (0, L): 69 69 if ((longitude[k] >= 0.) & (longitude[k] <= 90.)): # 4eme quadrant … … 90 90 ix[i] = 0 91 91 else: 92 ix[i] = math.ceil((x[i] - x0) / dx) - 1 # associates each x va kue to a cell number92 ix[i] = math.ceil((x[i] - x0) / dx) - 1 # associates each x value to a cell number 93 93 i = 0 94 94 iy = np.zeros([L], int) … … 101 101 # calculation of distances between point M(x,y) and 4 grid points of its belonging cell # 102 102 ######################################################################################### 103 close_point = np.zeros([L, 2], int) 103 close_point = np.zeros([L, 2], int) # for x- and y-axis, 2D-array of closest grid cell of new cartesian grid to observation 104 104 k = 0 105 105 for k in range (0, L): # boucle sur les points M (x et y) … … 108 108 d3 = sqrt(((x[k] - xvec[ix[k]]) ** 2) + ((y[k] - yvec[iy[k] + 1]) ** 2)) 109 109 d4 = sqrt(((x[k] - xvec[ix[k] + 1]) ** 2) + ((y[k] - yvec[iy[k] + 1]) ** 2)) 110 d_vec = np.array([d1, d2, d3, d4]) 110 d_vec = np.array([d1, d2, d3, d4]) # vector of distances between observation and each closer grid_cell of new cartesian grid 111 111 ind = np.where(d_vec == min(d_vec)) # finds the smallest distance between the 4 points of the grid 112 112 point_vec = np.array([(ix[k], iy[k]), (ix[k] + 1, iy[k]), (ix[k], iy[k] + 1), (ix[k] + 1, iy[k] + 1)])
Note: See TracChangeset
for help on using the changeset viewer.
