I’m translating some old IDL code into python. I’m very new and hoped (am still hoping) that by jumping into some IDL code I could maybe learn the language faster. I’ve started with what seemed to be the “easiest” function to bring over. This function takes a given “psd” and creates a simulated map from it. Not exactly sure what ‘psd’ stands for but we’re in an astronomical context so its probably safe to assume it means power spectral density… for our purposes this is a 2D array of whatever shape. Before taking some code that generates a fake atmosphere psd, I was using smaller 2×3, 3×4, & 4×4 arrays. Things were going… not exactly well but I was making it further in the code before things broke than I am now. I’ve got pretty limited experience with both IDL and numpy (astropy usually handled a lot of the numpy stuff for me so I never had to directly code new stuff), and so I’m really struggling to figure out just WHAT is wrong. The problem with fixing this arises because I have NO CLUE what parts of the IDL code need to stay in python.

So, the following is IDL code and my attempt to translate said IDL code into python.

IDL:

;we've allowed the map_ft too much freedom, so we need to make
;sure it will produce a real map

;obviously the dc component is real
map_ft[0,0] = sign(real_part(map_ft[0,0])) * abs(map_ft[0,0])

;if map_dim1 is divisible by 2 then we need to reverse
;along the extra column at the highest frequency
if (map_dim1 mod 2L) eq 0L then begin
  map_ft[map_dim1/2L,0L] = sign(real_part(map_ft[map_dim1/2L,0L])) * $
    abs(map_ft[map_dim1/2L,0L])
  map_ft[map_dim1/2L,1L:(map_dim2/2L - 1L + (map_dim2 mod 2L))] = $
    conj(reverse(reform(map_ft[map_dim1/2L,map_dim2/2L+1L:*])))
endif

Python:

# we've allowed the map_ft too much freedom, so we need to make
# sure it will produce a real map

# obviously the dc component is real (no clue what dc means btw)

map_ft[0, 0] = np.sign(map_ft[0, 0].real) * np.abs(map_ft[0, 0])

print(row_len, ' // 2 + (', row_len % 2, ') = ', (row_len // 2 + (row_len % 2)))

#if col_len is divisible by 2 then we need to reverse
#along the extra column at the highest frequency
if (col_len % 2) == 0: # working on columns 
    mid_col = col_len // 2
    map_ft[0, mid_col] = np.sign(map_ft[0, mid_col].real) * np.abs(map_ft[0, mid_col])
    

    end_row = (row_len // 2 + row_len % 2)
    map_ft[1 :row_len // 2 - 1 + (row_len % 2), mid_col] = 
        np.conj(np.flip(np.squeeze(map_ft[row_len // 2 + 1 :, mid_col])))

Note: map_dim2, mapdim1 = row_len, col_len = psd.shape & map_ft is a complex array of equal shape to psd (in this particular case psd is an array of size (4096, 2048)

the error occurs at the np.conj(…) line and reads: ValueError: could not broadcast input array from shape (2047,) into shape (2046,).

Finally we can get to what exactly I need help with: I’m wondering how necessary it is to translate the IDL code verbatim or if there’s more pythonic or cleaner ways of writing this code. I’m hesitant to try and simplify things just myself though because with my limited knowledge this would effectively mean abandoning the original code ad trying to understand the much harder to learn mathematical concepts that this code was founded on (end goal is to get an old atmospheric filtering code up and running, so lots of fft’s and scary math I haven’t learned yet. Also before anyone asks, I don’t have an IDL license (yet at least) and unfortunately can’t check what the code does when being run in IDL. For now just pretend its impossible and hopefully that changes.

I would appreciate help in making this code more pythonic because I learn quite quickly from working backwards and any help translating this code would be greatly appreciated.