Rebecca M. - 9 months ago 20

Python Question

I am trying to copy an extensive code from IDL into Python. One issue I believe I am having is with the definition of the function gridgen. Gridgen is a function used to generate a vertical grid with equal log-spaced grid points, where: zmin = coordinate at top of grid; zmax = coordinate at bottom of grid; Nlvls = desired number of levels in grid and z = output grid.

The IDL code is:

`FUNCTION GRIDGEN, zmin, zmax, Nlvls`

dlnz = (ALOG(zmax)-ALOG(zmin))/(Nlvls-1) ;divisions in log space

z = FLTARR(Nlvls) ;array of Nlvls points for logarithm to base 10

z[*] = 0. ;initialize grid to zero

z[Nlvls-1] = ALOG(zmax) ;assign the maximum value in log spacing

FOR i=Nlvls-2, 0, -1 DO z[i] = z[i+1] - dlnz ;generate log spacing values

z = EXP(z) ;convert from log(z) to actual values

RETURN, z

END

How I translated that into Python is:

`def gridgen100(zmin, zmax, Nlvls):`

dlnz = ((np.log(zmax) - np.log(zmin))/(Nlvls - 1)) # divisions in log space

z = np.zeros(Nlvls, dtype=float) # array of Nlvls points for logarithm to base 10

z[Nlvls-1] = np.log(zmax) # assign the maximum value in log spacing

for i in np.arange(Nlvls-2, Nlvls-101, -1): # NOT CORRECT; correct is: for i in [Nlvls-2, 0, -1]:

z[i] = z[i +1] - dlnz # generate log spacing values

#z = np.exp(np.array(z)) # convert from log(z) to actual values [MUST DO OUTSIDE DEF]

return z

The issues are:

- Because of how I created the for loop using np.arange, I have to define a separate gridgen function when I have different Nlvls (for example, sometimes I have 100 Nlvls, sometimes I have 25).
- I can't convert from log(z) to actual values within the function, I have to do it outside the definition.

I don't currently have access to IDL, so I'm unable to troubleshoot by comparing the IDL output to the Python output.

I am self-taught in Python and a beginner, but I appreciate any help or advice anyone can offer.

Answer

IIUC, your IDL for loop translates to

```
for i in range(Nlvls-2, -1, -1):
```

i.e. start at Nlvls-2 and drop by 1 until you reach 0. This gives me

```
def gridgen(zmin, zmax, Nlvls):
dlnz = (np.log(zmax) - np.log(zmin))/(Nlvls-1)
z = np.zeros(Nlvls, dtype=float)
z[Nlvls-1] = np.log(zmax)
for i in range(Nlvls-2, -1, -1):
z[i] = z[i+1] - dlnz
z = np.exp(z)
return z
```

and

```
>>> gridgen(2, 8, 10)
array([ 2. , 2.33305808, 2.72158 , 3.1748021 , 3.70349885,
4.32023896, 5.0396842 , 5.87893797, 6.85795186, 8. ])
```

But there's already a numpy function, `np.logspace`

, which does this log spacing for you, so if I'm right about what you're after, you could obtain the same result using it instead:

```
>>> np.logspace(np.log10(2), np.log10(8), 10)
array([ 2. , 2.33305808, 2.72158 , 3.1748021 , 3.70349885,
4.32023896, 5.0396842 , 5.87893797, 6.85795186, 8. ])
```

(For point #2, you can obviously remove the `z = np.exp(z)`

line if you don't want to return to the original space.)