How can I increase the speed of the following code with four for loops using other approaches like structure ?

D=dir('*.nc');

ERA_5=[];

for i=1:length(D)

D2=D(i);

name=D2.name;

Lon = ncread(name,"longitude");

Lat = ncread(name,"latitude");

Timee = ncread(name,"time");

T = datetime(1900,1,1) + hours(Timee);

T.Format = 'yyyy-MM-dd';

Wind_u= ncread(name,"u10");%with size =length(loon)* length(laat)* length(day)

Wind_v = ncread(name,"v10");%with size =length(loon)* length(laat)* length(day)

soil_temp = ncread(name,"stl1");%with size =length(loon)* length(laat)* length(day)

T_Precipitation = ncread(name,"tp");%with size =length(loon)* length(laat)* length(day)

for loon=1:length(Lon)

for laat=1:length(Lat)

for day=1:length(T)

T2=T(day);

[ y , m , d ] = ymd( T2);

ERA_5=[ERA_5; Lon(loon) Lat(laat) y m d Wind_u(loon, laat, day) Wind_v(loon, laat, day) soil_temp(loon, laat, day) T_Precipitation(loon, laat, day)];

end

end

end

end

[Harsh Saxena]

Try to use the code option in the editor while adding code to your question. Also, it will be better to if you attach the required files to run the code.

Answer 1

[Harsh Saxena]

Hi Mina,

The given code can be optimized in the following way:

1. If you can, prelocate the ERA_5 array, this will eliminate resizing of this array. If the size is not fixed, make it equal to a upper bound. For example, if length(D) is N, and the maximum dimensions of Lon, Lat, and T are Mlon, Mlat, and Mtime respectively, you can preallocate ERA_5 using ERA_5 = zeros(N*Mlon*Mlat*Mtime, 9);.
2. Use indexing instead of nested loops. You can create index variables idx_lon, idx_lat, and idx_day, which will increment inside the loops and determine the position in the ERA_5 array for each value. Initialize these variables before the loops with 1, and increment them at the end of each iteration (idx_day = idx_day + 1;, etc.).
3. Use the colon operator to read variables: Instead of using the ncread function inside the loops, you can use the colon operator : to read the entire variable in one go. This will avoid repeated function calls and improve performance. For example, replace Wind_u(loon, laat, day) with Wind_u(:, :, day).

So, the optimized version of your code can be like this:

D = dir('*.nc');
N = length(D);
Mlon = 0;
Mlat = 0;
Mtime = 0;

for i = 1:N
    D2 = D(i);
    name = D2.name;
    Lon = ncread(name, 'longitude');
    Lat = ncread(name, 'latitude');
    Timee = ncread(name, 'time');
    T = datetime(1900, 1, 1) + hours(Timee);
    T.Format = 'yyyy-MM-dd';
    
    if length(Lon) > Mlon
        Mlon = length(Lon);
    end
    
    if length(Lat) > Mlat
        Mlat = length(Lat);
    end
    
    if length(T) > Mtime
        Mtime = length(T);
    end
end

ERA_5 = zeros(N * Mlon * Mlat * Mtime, 9);
idx_lon = 1;
idx_lat = 1;
idx_day = 1;

for i = 1:N
    D2 = D(i);
    name = D2.name;
    Lon = ncread(name, 'longitude');
    Lat = ncread(name, 'latitude');
    Timee = ncread(name, 'time');
    T = datetime(1900, 1, 1) + hours(Timee);
    T.Format = 'yyyy-MM-dd';
    Wind_u = ncread(name, 'u10');
    Wind_v = ncread(name, 'v10');
    soil_temp = ncread(name, 'stl1');
    T_Precipitation = ncread(name, 'tp');
    
    for day = 1:length(T)
        T2 = T(day);
        [y, m, d] = ymd(T2);
        
        ERA_5(idx_lon:idx_lon+length(Lon)-1, 1) = Lon;
        ERA_5(idx_lat:idx_lat+length(Lat)-1, 2) = Lat;
        ERA_5(idx_day:idx_day+length(T)-1, 3) = y;
        ERA_5(idx_day:idx_day+length(T)-1, 4) = m;
        ERA_5(idx_day:idx_day+length(T)-1, 5) = d;
        ERA_5(idx_day:idx_day+length(T)-1, 6) = Wind_u(:, :, day);
        ERA_5(idx_day:idx_day+length(T)-1, 7) = Wind_v(:, :, day);
        ERA_5(idx_day:idx_day+length(T)-1, 8) = soil_temp(:, :, day);
        ERA_5(idx_day:idx_day+length(T)-1, 9) = T_Precipitation(:, :, day);
        
        idx_day = idx_day + length(T);
    end
    
    idx_lon = idx_lon + length(Lon);
    idx_lat = idx_lat + length(Lat);
end

% Trim the excess rows from the preallocated ERA_5 array
ERA_5 = ERA_5(1:idx_day-1, :);

Note that the code assumes that the dimensions of Lon, Lat, and T remain constant across all files. If this is not the case, you may need to modify the code accordingly.

Also note that this is an intuition based code and some errors might occur which will be trivial to resolve. You may focus on the intuition rather than the actual code.

Hope this helps!

[Mina Mino]

Hi Harsh Saxena,

Thank you so much for your consideration and answer!

Best,

Mina