how can I reduce the computational time of a single for loop in a Matlab code

Question

Mina Mino 2023-4-8

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1943644-how-can-i-reduce-the-computational-time-of-a-single-for-loop-in-a-matlab-code

编辑：埃博拉酱 2023-4-9

I have a problem with the following for. it is too time consuming. could anyone help me with it?thank in advance for your consideration and time.

Let me give you some more details about the aim of the code. I have observations with the number of 'd' that measured in 'day-of-year'. so each observations occured in a specific 'day-of-year'.

I want to collocate each of these observations (10 million data) with 'VWC', 'SM_R' that measured in coordinates (Lat, Lon) and (LAT_R, LON_R) in days called 'DOY' and 'DOY-R' in code, respectivelly. see the associate part of code:

%%%%%Matlab Code%%%%%%

Lon_d=Lon(DOY==day_of_year(d));

Lat_d=Lat(DOY==day_of_year(d));

VWC_SMAP_d=VWC(DOY==day_of_year(d));

LAT_R_d=LAT_R(DOY_R==day_of_year(d));

LON_R_d= LON_R(DOY_R==day_of_year(d));

SM_R_9_d=SM_R(DOY_R==day_of_year(d));

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

to achieve this goal, each observation with coordinate 'lon_cy(d), lat_cy(d)' was collocated with observations limited to 'step (18/110 degree) *step (18/110 degree)' grids with the center of it (lon_cy(d), lat_cy(d)):

%%%%%%%Matlab Code%%%%%%%%%%%%%%%%%%%%

lon_smap_9_index= find(LON_R_d <=lon_cy(d)+step & LON_R_d>=lon_cy(d)-step);

lat_smap_9_index= find(LAT_R_d <=lat_cy(d)+step & LAT_R_d>=lat_cy(d)-step);

intsec_smap_9=intersect(lon_smap_9_index,lat_smap_9_index);

VWC_SMAP_g_9=VWC_R_9_d(intsec_smap_9)

lon_smap_36_index= find(Lat_d <=lat_cy (d)+step & Lat_d>=lat_cy(d)-step);

lat_smap_36_index= find(Lon_d <=lon_cy(d)+step & Lon_d>=lon_cy(d)-step);

intsec_smap_36=intersect(lon_smap_36_index,lat_smap_36_index);

SM_g_36=SM_d(intsec_smap_36);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

two additional ancillary datasets were aimed to colloated with the d-dimension observations for each day-of-year:

here I just need to find the nearest ancillary datasets to each observation. that is why I used 'knnsearch' function.

%%%%%%%%%%%%%%%%%%%Matlab%%%%%%%%%%%%%%

point_1=[lat_cy(d) lon_cy(d)];

ptCloud_1=[latitu longi];

idx_land = knnsearch(ptCloud_1,point_1);

ptCloud_2=[clay(:,1) clay(:,2)];

idx_clay = knnsearch(ptCloud_2,point_1);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

and finally I will store each collocated observations in the following way and save it

%%%%%%%%%%%%%%%%%%%%Matlab%%%%%%%%%%%%%%%

FT=[day_of_year(d) lon_cy(d) lat_cy(d) inc_cy(d) ref_cy(d) mean(VWC_SMAP_g_9) mean(roughness_g_9) mean(SM_g_9) mean(VWC_SMAP_g_36) mean(roughness_g_36) mean(TM_d_36) mean(SM_g_36) clay(idx_clay,3) land_cover(idx_land)];

8 个评论
显示 6更早的评论隐藏 6更早的评论

Mina Mino 2023-4-8

编辑：Walter Roberson 2023-4-8

在 MATLAB Online 中打开

thanks for your reply. I changed it in the following form. But is does not improve significantly. could you please give more tips? Thanks in advance.

%%
FT=zeros(1,14);step=18/110;
parfor  d=8350381:10350380 
    
    
    Lon_d=Lon(DOY==day_of_year(d));
    Lat_d=Lat(DOY==day_of_year(d));
    VWC_SMAP_d=VWC(DOY==day_of_year(d));
    roughness_d=roughness(DOY==day_of_year(d));
    TM_d=TM(DOY==day_of_year(d));
    SM_d=SM_36(DOY==day_of_year(d));
    
    LAT_R_d=LAT_R(DOY_R==day_of_year(d));
    LON_R_d= LON_R(DOY_R==day_of_year(d));
    VWC_R_9_d=VWC_R(DOY_R==day_of_year(d));  
    SM_R_9_d=SM_R(DOY_R==day_of_year(d));
    ROUGH_R_9_d=ROUGH_R(DOY_R==day_of_year(d));
    
    if (~isempty(ROUGH_R_9_d) && ~isempty(SM_d))
        
        lon_smap_9_index= find(LON_R_d <=lon_cy(d)+step & LON_R_d>=lon_cy(d)-step); 
        lat_smap_9_index= find(LAT_R_d <=lat_cy(d)+step & LAT_R_d>=lat_cy(d)-step);
        intsec_smap_9=intersect(lon_smap_9_index,lat_smap_9_index);
        
        
        lon_smap_36_index= find(Lat_d <=lat_cy (d)+step & Lat_d>=lat_cy(d)-step);
        lat_smap_36_index= find(Lon_d <=lon_cy(d)+step & Lon_d>=lon_cy(d)-step);
        intsec_smap_36=intersect(lon_smap_36_index,lat_smap_36_index); 
        
        
        VWC_SMAP_g_9=VWC_R_9_d(intsec_smap_9);roughness_g_9=ROUGH_R_9_d(intsec_smap_9);SM_g_9=SM_R_9_d(intsec_smap_9);
        
        VWC_SMAP_g_9(isnan(VWC_SMAP_g_9))=[]; 
        roughness_g_9(isnan(roughness_g_9))=[]; 
        SM_g_9(isnan(SM_g_9))=[]; 
        
        
        VWC_SMAP_g_36=VWC_SMAP_d(intsec_smap_36);roughness_g_36=roughness_d(intsec_smap_36);TM_d_36=TM_d(intsec_smap_36) ;SM_g_36=SM_d(intsec_smap_36);
        
        VWC_SMAP_g_36(isnan(VWC_SMAP_g_36))=[]; 
        roughness_g_36(isnan(roughness_g_36))=[];
        TM_d_36(isnan(TM_d_36))=[];
        SM_g_36(isnan(SM_g_36))=[]; 
        
        
        point_1=[lat_cy(d) lon_cy(d)];
        ptCloud_1=[latitu longi];
        idx_land  = knnsearch(ptCloud_1,point_1);
        
        ptCloud_2=[clay(:,1) clay(:,2)];
        idx_clay  = knnsearch(ptCloud_2,point_1);
        
        
        if (~isempty(intsec_smap_36)&& ~isempty(intsec_smap_9)) 
            
            FT=[day_of_year(d) lon_cy(d) lat_cy(d) inc_cy(d) ref_cy(d) mean(VWC_SMAP_g_9) mean(roughness_g_9) mean(SM_g_9) mean(VWC_SMAP_g_36) mean(roughness_g_36) mean(TM_d_36) mean(SM_g_36) clay(idx_clay,3) land_cover(idx_land)];
            S=FT;
            parsave(sprintf('output%d.mat', d),S );
        end
    end
end

Mina Mino 2023-4-8

Let me give you some more details about the aim of the code. I have observations with the number of 'd' that measured in 'day-of-year'. so each observations occured in a specific 'day-of-year'.

I want to collocate each of these observations (10 million data) with 'VWC', 'SM_R' that measured in coordinates (Lat, Lon) and (LAT_R, LON_R) in days called 'DOY' and 'DOY-R' in code, respectivelly. see the associate part of code:

%%%%%Matlab Code%%%%%%

Lon_d=Lon(DOY==day_of_year(d));

Lat_d=Lat(DOY==day_of_year(d));

VWC_SMAP_d=VWC(DOY==day_of_year(d));

LAT_R_d=LAT_R(DOY_R==day_of_year(d));

LON_R_d= LON_R(DOY_R==day_of_year(d));

SM_R_9_d=SM_R(DOY_R==day_of_year(d));

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

to achieve this goal, each observation with coordinate 'lon_cy(d), lat_cy(d)' was collocated with observations limited to 'step (18/110 degree) *step (18/110 degree)' grids with the center of it (lon_cy(d), lat_cy(d)):

%%%%%%%Matlab Code%%%%%%%%%%%%%%%%%%%%

lon_smap_9_index= find(LON_R_d <=lon_cy(d)+step & LON_R_d>=lon_cy(d)-step);

lat_smap_9_index= find(LAT_R_d <=lat_cy(d)+step & LAT_R_d>=lat_cy(d)-step);

intsec_smap_9=intersect(lon_smap_9_index,lat_smap_9_index);

VWC_SMAP_g_9=VWC_R_9_d(intsec_smap_9)

lon_smap_36_index= find(Lat_d <=lat_cy (d)+step & Lat_d>=lat_cy(d)-step);

lat_smap_36_index= find(Lon_d <=lon_cy(d)+step & Lon_d>=lon_cy(d)-step);

intsec_smap_36=intersect(lon_smap_36_index,lat_smap_36_index);

SM_g_36=SM_d(intsec_smap_36);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

two additional ancillary datasets were aimed to colloated with the d-dimension observations for each day-of-year:

here I just need to find the nearest ancillary datasets to each observation. that is why I used 'knnsearch' function.

%%%%%%%%%%%%%%%%%%%Matlab%%%%%%%%%%%%%%

point_1=[lat_cy(d) lon_cy(d)];

ptCloud_1=[latitu longi];

idx_land = knnsearch(ptCloud_1,point_1);

ptCloud_2=[clay(:,1) clay(:,2)];

idx_clay = knnsearch(ptCloud_2,point_1);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

and finally I will store each collocated observations in the following way and save it

%%%%%%%%%%%%%%%%%%%%Matlab%%%%%%%%%%%%%%%

FT=[day_of_year(d) lon_cy(d) lat_cy(d) inc_cy(d) ref_cy(d) mean(VWC_SMAP_g_9) mean(roughness_g_9) mean(SM_g_9) mean(VWC_SMAP_g_36) mean(roughness_g_36) mean(TM_d_36) mean(SM_g_36) clay(idx_clay,3) land_cover(idx_land)];

Walter Roberson 2023-4-8

在 MATLAB Online 中打开

    Lon_d=Lon(DOY==day_of_year(d));
    Lat_d=Lat(DOY==day_of_year(d));
    VWC_SMAP_d=VWC(DOY==day_of_year(d));
    roughness_d=roughness(DOY==day_of_year(d));
    TM_d=TM(DOY==day_of_year(d));
    SM_d=SM_36(DOY==day_of_year(d));
    
    LAT_R_d=LAT_R(DOY_R==day_of_year(d));
    LON_R_d= LON_R(DOY_R==day_of_year(d));
    VWC_R_9_d=VWC_R(DOY_R==day_of_year(d));  
    SM_R_9_d=SM_R(DOY_R==day_of_year(d));
    ROUGH_R_9_d=ROUGH_R(DOY_R==day_of_year(d));

You can pre-calculate those for each different value 1 to 365 (or 366 if appropriate), and index them by day_of_year(d), instead of calculating them every time.

Mina Mino 2023-4-8

编辑：Mina Mino 2023-4-8

thanks for your reply. I rewrite the code in the following form, the elapsed time decreased from 2.3 to 0.3 second for each day-of-year. However, it also takes a days and even months for 12 million or 20 milion observations. I think I have to change the way of coding. several days ago one of our friends here recommend a way of thinking for another code which reduce the computational time from months to a few second. Do you have any idea?

%%

FT=zeros(1,14);step=18/110;

SMAP_36=[DOY Lon Lat VWC roughness TM SM_36 ];

SMAP_9=[DOY_R LON_R LAT_R VWC_R ROUGH_R SM_R];

parfor d=8350381:10350380

tic

SMAP_36_d=SMAP_36(SMAP_36(:,1)==day_of_year(d),:);

SMAP_9_d=SMAP_9(SMAP_9(:,1)==day_of_year(d),:);

if (~isempty(SMAP_36_d) && ~isempty(SMAP_9_d))

lon_smap_9_index= find(SMAP_36_d(:,2) <=lon_cy(d)+step & SMAP_36_d(:,2)>=lon_cy(d)-step);

lat_smap_9_index= find(SMAP_36_d(:,3) <=lat_cy(d)+step & SMAP_36_d(:,3)>=lat_cy(d)-step);

intsec_smap_9=intersect(lon_smap_9_index,lat_smap_9_index);

lon_smap_36_index= find(SMAP_9_d(:,3) <=lat_cy (d)+step & SMAP_9_d(:,3)>=lat_cy(d)-step);

lat_smap_36_index= find(SMAP_9_d(:,2) <=lon_cy(d)+step & SMAP_9_d(:,2)>=lon_cy(d)-step);

intsec_smap_36=intersect(lon_smap_36_index,lat_smap_36_index);

if (~isempty(intsec_smap_36)&& ~isempty(intsec_smap_9))

VWC_SMAP_g_9=SMAP_9(intsec_smap_9,4);roughness_g_9=SMAP_9(intsec_smap_9,5);SM_g_9=SMAP_9(intsec_smap_9,6);

VWC_SMAP_g_9(isnan(VWC_SMAP_g_9))=[];

roughness_g_9(isnan(roughness_g_9))=[];

SM_g_9(isnan(SM_g_9))=[];

VWC_SMAP_g_36=SMAP_36(intsec_smap_36,4);roughness_g_36=SMAP_36(intsec_smap_36,5);TM_d_36=SMAP_36(intsec_smap_36,6) ;SM_g_36=SMAP_36(intsec_smap_36,7);

VWC_SMAP_g_36(isnan(VWC_SMAP_g_36))=[];

roughness_g_36(isnan(roughness_g_36))=[];

TM_d_36(isnan(TM_d_36))=[];

SM_g_36(isnan(SM_g_36))=[];

point_1=[lat_cy(d) lon_cy(d)];

ptCloud_1=[latitu longi];

idx_land = knnsearch(ptCloud_1,point_1);

ptCloud_2=[clay(:,1) clay(:,2)];

idx_clay = knnsearch(ptCloud_2,point_1);

FT=[day_of_year(d) lon_cy(d) lat_cy(d) inc_cy(d) ref_cy(d) mean(VWC_SMAP_g_9) mean(roughness_g_9) mean(SM_g_9) mean(VWC_SMAP_g_36) mean(roughness_g_36) mean(TM_d_36) mean(SM_g_36) clay(idx_clay,3) land_cover(idx_land)];

S=FT;

% parsave(sprintf('output%d.mat', d),S );

end

toc

end

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

埃博拉酱 2023-4-9

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1943644-how-can-i-reduce-the-computational-time-of-a-single-for-loop-in-a-matlab-code#answer_1212184

编辑：埃博拉酱 2023-4-9

It seems that you are calling knnsearch twice in each loop to search for a single point in the point cloud. However, knnsearch can actually search for multiple points at once. It is recommended that you try to move knnsearch out of the loop and use the vectorized syntax of knnsearch to calculate all idx_land and idx_clay at once.

Based on my experience, I feel that the bottleneck of your code's speed is this knnsearch. If it's not, please reply to me again.

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

how can I reduce the computational time of a single for loop in a Matlab code

8 个评论
显示 6更早的评论隐藏 6更早的评论

回答（1 个）

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

how can I reduce the computational time of a single for loop in a Matlab code

8 个评论 显示 6更早的评论隐藏 6更早的评论

回答（1 个）

0 个评论 显示 -2更早的评论隐藏 -2更早的评论

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

8 个评论
显示 6更早的评论隐藏 6更早的评论

0 个评论
显示 -2更早的评论隐藏 -2更早的评论