function [ind_obs,ind_bpj,ind_val] = gen_warp_mat(h,w,H,ih)   

% from myWarp.m
[x,y] = meshgrid(1:w,ih);
pix   = [x(:)'; y(:)'];
hPixels = [ pix; ones(1,w)];
hScene  = H*hPixels; % back prj observed pixel grid to ref grid (H is in fact the camera motion)
loc_x=(hScene(1,:)./(hScene(3,:)));
loc_x_floor = floor(loc_x);
loc_y=(hScene(2,:)./(hScene(3,:)));
loc_y_floor = floor(loc_y); 
% end from myWarp.m

dif1 = loc_x - loc_x_floor;
dif2 = loc_y - loc_y_floor;
weight = [(1-dif1).*(1-dif2), (1-dif1).*dif2, dif1.*(1-dif2), dif1.*dif2];

temp_loc = (loc_x_floor-1)*h+loc_y_floor;
temp_bckprjed_ind = [temp_loc, temp_loc+1, temp_loc+h, temp_loc+h+1];
temp_obsvd_ind = [ih:h:h*w, ih:h:h*w, ih:h:h*w, ih:h:h*w];
temp_value_ind = weight;

temp_warped_x = [loc_x_floor, loc_x_floor, loc_x_floor, loc_x_floor];
temp_warped_y = [loc_y_floor, loc_y_floor, loc_y_floor, loc_y_floor];
final_ind = find( temp_warped_x>0 & temp_warped_x<w & temp_warped_y>0 & temp_warped_y<h );

ind_obs = temp_obsvd_ind(final_ind);
ind_bpj = temp_bckprjed_ind(final_ind);
ind_val = temp_value_ind(final_ind);

% mtx_warp = sparse(ind_obs,ind_bpj,ind_val,h*w,h*w);
