Ndtyjky

I have this code that I wrote in MATLAB to store matrices. I was using cells array, but in Python I don't know how to do this.

Would someone know how can I do this?

The MATLAB code is:

S =[0.5 0.7 0.9 1.1]; % distância entre tx e rx[m]

d = 0.2*ones(1,10); 

h = [ 0 0.1 0.2 0.3 0.4]



n_cam = length(d); % numero de camadas

n_alt = length(h); % numero de alturas

n_S = length(S); % numero de receptores

z = zeros(n_cam,n_alt); % profundidade

Rv_h = zeros(n_S,n_alt);

Rv_v = zeros(n_S,n_alt);

Rv = zeros(n_cam,n_alt);

Rh = zeros(n_cam,n_alt);

S_Rv = cell(1,n_S);

S_Rh = cell(1,n_S);

sigma = 0.3*ones(1,n_cam);

sigmaah = zeros(n_S,n_alt);





for i = 1:n_S

    for  j = 1:n_alt

        for k = 1:n_cam

            z(k,j)= (sum(d(1:k))+h(j))/S(i);

            Rv(k,j) = 1/((4*z(k,j)^2+1)^0.5);

            Rh(k,j) = ((4*z(k,j)^2+1)^0.5)-2*z(k,j);

        end 

        Rv_h(i,j) = 1/((4*(h(j)/S(i))^2+1)^0.5);

        Rh_h(i,j)=((4*(h(j)/S(i))^2+1)^0.5)-2*(h(j)/S(i));

    end

    S_Rv(:,i) = {Rv}; % z para cada camada em cada altura, para cada S

    S_Rh(:,i) = {Rh};

end



for i = 1:n_S

    for  j = 1:n_alt

        Rv = cell2mat(S_Rv(1,i));

        Rh = cell2mat(S_Rh(1,i));

        sigma_ah = sigma(1)*(Rh_h(i,j)-Rh(1,j));

        sigma_av = sigma(1)*(Rv_h(i,j)-Rv(1,j));

        for k = 2:(n_cam-1)

           sigma_ah_ant = sigma_ah;

           sigma_av_ant = sigma_av;

           sigma_ah = sigma_ah_ant + sigma(k)*(Rh(k-1,j)-Rh(k,j));

           sigma_av = sigma_av_ant + sigma(k)*(Rv(k-1,j)-Rv(k,j));

        end

        sigmaah (i,j)  = sigma_ah + sigma(end)*Rh(n_cam-1,j)

        sigmaav (i,j)  = sigma_av + sigma(end)*Rv(n_cam-1,j)

    end

end

I was thinking that in Python I could do something like:

n_S = 4

n_alt = 9

n_cam = 6

Rv =

for i in range(1,n_S):

    for j in range(1,n_alt):

        for k in range(1,n_cam):

            z[k][j]= (sum(d[0:k])+h[j])/S[i]

            Rv[i][j][k] = 1/((4*z[k,j]**2+1)**0.5)

But it is not working, the error message I get is "list index out of range".

I would suggest taking use of the broadcasting feature

The loops can then be replaced to clarify the code :

from numpy import array, ones

S =array([0.5, 0.7, 0.9, 1.1])

d = 0.2*ones((10)); 

h = array([ 0, 0.1, 0.2, 0.3, 0.4])



z = ((d.cumsum()[:, None] + h).ravel() / S[:, None]).reshape((S.size, d.size, h.size))

Rv = 1 / (4 * z ** 2 + 1)** .5

# ...etc

print(Rv[-1])

Which outputs:

    [[0.93979342 0.87789557 0.80873608 0.73994007 0.67572463]

 [0.80873608 0.73994007 0.67572463 0.61782155 0.56652882]

 [0.67572463 0.61782155 0.56652882 0.52145001 0.48191875]

 [0.56652882 0.52145001 0.48191875 0.4472136  0.41665471]

 [0.48191875 0.4472136  0.41665471 0.38963999 0.36565237]

 [0.41665471 0.38963999 0.36565237 0.34425465 0.32507977]

 [0.36565237 0.34425465 0.32507977 0.30782029 0.29221854]

 [0.32507977 0.30782029 0.29221854 0.27805808 0.26515648]

 [0.29221854 0.27805808 0.26515648 0.25335939 0.24253563]

 [0.26515648 0.25335939 0.24253563 0.23257321 0.22337616]]

Which overlaps with the computation in octave/matlab:

Rv =



   0.93979   0.87790   0.80874   0.73994   0.67572

   0.80874   0.73994   0.67572   0.61782   0.56653

   0.67572   0.61782   0.56653   0.52145   0.48192

   0.56653   0.52145   0.48192   0.44721   0.41665

   0.48192   0.44721   0.41665   0.38964   0.36565

   0.41665   0.38964   0.36565   0.34425   0.32508

   0.36565   0.34425   0.32508   0.30782   0.29222

   0.32508   0.30782   0.29222   0.27806   0.26516

   0.29222   0.27806   0.26516   0.25336   0.24254

   0.26516   0.25336   0.24254   0.23257   0.22338   

This reduces the pyramid of doom and is probably faster than the for-loops due to numpy magic.
Edit: formatting
Edit2: gave check

You need to define z and Rv to be 2D and 3D arrays of known size before you do this.

Note also that python (and numpy) arrays are zero-based, so just use range directly.

Not tested:

import numpy as np



d = 0.2*np.arange(10)   ### not sure if this is what you meant

h = np.array([ 0 0.1 0.2 0.3 0.4])



n_S = 4

n_alt = 9

n_cam = 6

Rv = np.zeros((n_S, n_alt, n_cam))

z = np.zeros((n_cam, n_alt))

for i in range(n_S):

    for j in range(n_alt):

        for k in range(n_cam):

            z[k][j]= (sum(d[0:k])+h[j])/S[i]   ## not sure about 0:10

            Rv[i][j][k] = 1/((4*z[k,j]**2+1)**0.5)

However, as @GlobalTraveler points out, the most pythonic/numpyish way to do this is to take advantage of broadcasting and not use a loop at all:

Rv = 1/np.sqrt(4 * z**2 + 1)