Changeset 1160 for trunk/src/series/civ_series.m
- Timestamp:
- Jul 16, 2024, 9:57:54 AM (4 months ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/src/series/civ_series.m
r1158 r1160 952 952 % .DVDY 953 953 954 function [xtable,ytable,utable,vtable,ctable,F ,result_conv,errormsg] = civ (par_civ)954 function [xtable,ytable,utable,vtable,ctable,FF,result_conv,errormsg] = civ (par_civ) 955 955 956 956 %% prepare measurement grid … … 972 972 par_civ.Grid(:,2)=reshape(GridY,[],1);% increases with array index 973 973 % 974 %975 974 % minix=floor(par_civ.Dx/2)-0.5; 976 975 % maxix=minix+par_civ.Dx*floor((par_civ.ImageWidth-1)/par_civ.Dx); … … 1001 1000 vtable=shifty;%zeros(nbvec,1); 1002 1001 ctable=zeros(nbvec,1); 1003 F =zeros(nbvec,1);1002 FF=zeros(nbvec,1); 1004 1003 result_conv=[]; 1005 1004 errormsg=''; … … 1058 1057 iref=round(par_civ.Grid(ivec,1)+0.5);% xindex on the image A for the middle of the correlation box 1059 1058 jref=round(par_civ.ImageHeight-par_civ.Grid(ivec,2)+0.5);% j index for the middle of the correlation box in the image A 1060 F (ivec)=0;1059 FF(ivec)=0; 1061 1060 subrange1_x=iref-ibx2:iref+ibx2;% x indices defining the first subimage 1062 1061 subrange1_y=jref-iby2:jref+iby2;% y indices defining the first subimage … … 1078 1077 sizemask=sum(sum(mask1_crop))/(numel(subrange1_y)*numel(subrange1_x));%size of the masked part relative to the correlation sub-image 1079 1078 if sizemask > 1/2% eliminate point if more than half of the correlation box is masked 1080 F (ivec)=1; %1079 FF(ivec)=1; % 1081 1080 utable(ivec)=NaN; 1082 1081 vtable(ivec)=NaN; … … 1092 1091 end 1093 1092 %threshold on image minimum 1094 if F (ivec)~=11093 if FF(ivec)==0 1095 1094 if check_MinIma && (image1_mean < par_civ.MinIma || image2_mean < par_civ.MinIma) 1096 F (ivec)=1;1095 FF(ivec)=1; 1097 1096 %threshold on image maximum 1098 1097 elseif check_MaxIma && (image1_mean > par_civ.MaxIma || image2_mean > par_civ.MaxIma) 1099 F (ivec)=1;1100 end 1101 if F (ivec)==11098 FF(ivec)=1; 1099 end 1100 if FF(ivec)==1 1102 1101 utable(ivec)=NaN; 1103 1102 vtable(ivec)=NaN; … … 1127 1126 sum_square=sum(sum(image1_crop.*image1_crop)); 1128 1127 %reference: Oliver Pust, PIV: Direct Cross-Correlation 1128 %%%%%% correlation calculation 1129 1129 result_conv= conv2(image2_crop,flip(flip(image1_crop,2),1),'valid'); 1130 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 1130 1131 corrmax= max(max(result_conv)); 1131 1132 result_conv=(result_conv/corrmax)*255; %normalize, peak=always 255 … … 1138 1139 try 1139 1140 if par_civ.CorrSmooth==1 1140 [vector,F (ivec)] = SUBPIXGAUSS (result_conv,x,y);1141 [vector,FF(ivec)] = SUBPIXGAUSS (result_conv,x,y); 1141 1142 elseif par_civ.CorrSmooth==2 1142 [vector,F (ivec)] = SUBPIX2DGAUSS (result_conv,x,y);1143 [vector,FF(ivec)] = SUBPIX2DGAUSS (result_conv,x,y); 1143 1144 else 1144 [vector,F (ivec)] = quadr_fit(result_conv,x,y);1145 [vector,FF(ivec)] = quadr_fit(result_conv,x,y); 1145 1146 end 1146 1147 utable(ivec)=vector(1)*mesh+shiftx(ivec); … … 1154 1155 utable(ivec)=0; 1155 1156 vtable(ivec)=0; 1156 F (ivec)=1;1157 FF(ivec)=1; 1157 1158 end 1158 1159 ctable(ivec)=corrmax/sum_square;% correlation value 1159 1160 catch ME 1160 F (ivec)=1;1161 FF(ivec)=1; 1161 1162 disp(ME.message) 1162 1163 end 1163 1164 else 1164 F (ivec)=1;1165 FF(ivec)=1; 1165 1166 end 1166 1167 end … … 1188 1189 %http://urapiv.wordpress.com 1189 1190 peaky = y;peakx=x; 1190 if y < npy && y > 1 && x < npx-1 && x > 11191 f0 = log(result_conv(y,x));1192 f1 = log(result_conv(y-1,x));1193 f2 = log(result_conv(y+1,x));1194 peaky = peaky+ (f1-f2)/(2*f1-4*f0+2*f2);1195 f1 = log(result_conv(y,x-1));1196 f2 = log(result_conv(y,x+1));1197 peakx = peakx+ (f1-f2)/(2*f1-4*f0+2*f2);1198 else1199 F=1; % warning flag for vector truncated by the limited search box1200 end1191 % if y < npy && y > 1 && x < npx-1 && x > 1 1192 % f0 = log(result_conv(y,x)); 1193 % f1 = log(result_conv(y-1,x)); 1194 % f2 = log(result_conv(y+1,x)); 1195 % peaky = peaky+ (f1-f2)/(2*f1-4*f0+2*f2); 1196 % f1 = log(result_conv(y,x-1)); 1197 % f2 = log(result_conv(y,x+1)); 1198 % peakx = peakx+ (f1-f2)/(2*f1-4*f0+2*f2); 1199 % else 1200 % F=1; % warning flag for vector truncated by the limited search box 1201 % end 1201 1202 1202 1203 vector=[peakx-floor(npx/2)-1 peaky-floor(npy/2)-1];
Note: See TracChangeset
for help on using the changeset viewer.