diff --git a/src/ich/args_in.py b/src/ich/args_in.py index 5d0f46502ec87e907a664e10796bf0958a6613db..5c6aa30a11d7443a78a4c084859c8e7dd5767b21 100644 --- a/src/ich/args_in.py +++ b/src/ich/args_in.py @@ -31,11 +31,10 @@ np.seterr(divide='ignore', invalid='ignore') from plotting import density_scatter,plot_data,rmse -from algo import arc_sinc,arc_sinc_,cal_ +from algo import arc_sinc,arc_sinc_,cal_, dynamicWindow from rst_io import read_bin, write_tif, write_bin from utils import blockshaped,unblockshaped,unpad,spatial_intp_lin,split_sample ########################################################################## - def rvog_inverse(args): try: @@ -168,59 +167,12 @@ def rvog_inverse(args): else: - temp_cor = blockshaped(cor, args.window_size, args.window_size) - temp_lidar = blockshaped(lidar_ht_cal, args.window_size, args.window_size) - - s = np.zeros(np.shape(temp_cor)) - c = np.zeros(np.shape(temp_cor)) - ht_ = np.zeros(np.shape(temp_cor)) - rmse__ = np.zeros(np.shape(temp_cor)) - count = np.zeros(np.shape(temp_cor)) - - - parm_ = [0,0,0,0,0] - print("[1/3] Generating calibration parameters...") - for win in tqdm(range(np.shape(temp_cor)[0])): - # for win in (range(np.shape(temp_cor)[0])): - parm = cal_(temp_cor[win,:,:],temp_lidar[win,:,:],args.htl,args.htg) - - mask = temp_lidar[win,:,:].copy() - mask[~np.isnan(mask)] = 1 - if np.all(np.array(parm) == 0): - parm=parm_.copy() - if np.all(mask==0) or np.all(np.isnan(mask)): - # mask[np.shape(mask)[0]//2,np.shape(mask)[1]//2]=1 - np.fill_diagonal(mask, 1) - mask = np.flipud(mask) - np.fill_diagonal(mask, 1) - - - s[win,:,:] = np.full((args.window_size,args.window_size), parm[1])*mask - c[win,:,:] = np.full((args.window_size,args.window_size), parm[2])*mask - # s[win,:,:] = np.full((args.window_size,args.window_size), parm[1]) - # c[win,:,:] = np.full((args.window_size,args.window_size), parm[2]) - - rmse__[win,:,:] = np.full((args.window_size,args.window_size), parm[4]) - count[win,:,:] = np.full((args.window_size,args.window_size), - np.count_nonzero(~np.isnan(temp_lidar[win,:,:]))) - # print(np.count_nonzero(~np.isnan(temp_lidar[win,:,:])),parm[1],parm[2]) - gama = temp_cor[win,:,:] / parm[1] - ht_[win,:,:] = arc_sinc(gama, parm[2])*mask - parm_ = parm.copy() - - s = unblockshaped(s, rows,cols) - c = unblockshaped(c, rows,cols) - rmse__ = unblockshaped(rmse__, rows,cols) - ht_ = unblockshaped(ht_, rows,cols) - temp_cor = unblockshaped(temp_cor, rows,cols) - count = unblockshaped(count, rows,cols) - - + print("[1/3] Generating calibration parameters...") + s,c,rmse__,ht_, count = dynamicWindow(cor,lidar_ht_cal, args.window_size,args.htl,args.htg) ht_[ht_==0] = np.nan s[s==0]=np.nan c[c==0]=np.nan - rmse__[rmse__==0]=np.nan if rPad!=0 or cPad!=0: @@ -230,9 +182,12 @@ def rvog_inverse(args): ht_ = unpad(ht_, pad_width) rmse__ = unpad(rmse__, pad_width) count = unpad(count, pad_width) - temp_cor = unpad(temp_cor, pad_width) - + temp_cor = unpad(cor, pad_width) + print(lidar_ht_cal.shape,ht_.shape) + print(np.nanmean(lidar_ht_cal),np.nanstd(lidar_ht_cal)) + print(np.nanmean(ht_),np.nanstd(ht_)) + plot_data(lidar_ht_cal,ht_,args.htl,args.htg,'Lidar height(m)','InSAR inverted (m)',pltName)