%% Test locations - see...well see previous everything. % Mumbai 18.3N, 73.3E % Rajasthan %25N, 74E % Kerala %10N, 76.6E % Shimla %32.04N, 76.4E % Delhi %28.44N, 77.12E % Hyderabad %17.52N, 77.8E % Patna %26.4N, 84.8E % Bhubanswar %20.5N, 85.7E % Meghalaya %25.56N, 91.52E IMD_loc_lon_points=[29 32 42 41 44 47 75 78 101]; IMD_loc_lat_points=[48 75 15 103 89 45 81 57 77]; loc_opts={'Mumbai','Rajasthan','Kerala','Shimla','Delhi',... 'Hyderabad','Patna','Bhubaneswar','Meghalaya'}; %loop to pick out only those points. for loc=1:9 IMD_tps(loc,:)=squeeze(IMD_mapped_frdates(IMD_loc_lat_points(loc),... IMD_loc_lon_points(loc),:)); IMD_tps_all(loc,:)=squeeze(IMD_map_rorder(IMD_loc_lat_points(loc),... IMD_loc_lon_points(loc),:)); IMERG_tps(loc,:)=squeeze(IMERG_to_IMD_frdates(IMD_loc_lat_points(loc),... IMD_loc_lon_points(loc),:)); Fcst_tps(loc,:,:)=squeeze(fr12km_d1_ensall(IMD_loc_lat_points(loc),... IMD_loc_lon_points(loc),:,:)); UQM_tps(loc,:,:)=squeeze(IMD_QM_DGQM_ensall(IMD_loc_lat_points(loc),... IMD_loc_lon_points(loc),:,:)); EMOS_tps(loc,:,:)=squeeze(fr12km_EMOS_ensall(IMD_loc_lat_points(loc),... IMD_loc_lon_points(loc),:,:)); end % %add a final point, India mean. Need to reshape into latxlon matrix. % IMD_tps(10,:)=nanmean(reshape(IMD_mapped_frdates,... % length(IMD_lat_points)*length(IMD_lon_points),... % size(IMD_mapped_frdates,3)),1); % IMERG_tps(10,:)=nanmean(reshape(IMERG_to_IMD_frdates,... % length(IMD_lat_points)*length(IMD_lon_points),... % size(IMERG_to_IMD_frdates,3)),1); % Fcst_tps(10,:,:)=nanmean(reshape(fr12km_d1_ensall,... % length(IMD_lat_points)*length(IMD_lon_points),... % size(fr12km_d1_ensall,3),size(fr12km_d1_ensall,4)),1); % UQM_tps(10,:,:)=nanmean(reshape(IMD_QM_DGQM_ensall,... % length(IMD_lat_points)*length(IMD_lon_points),... % size(IMD_QM_DGQM_ensall,3),size(IMD_QM_DGQM_ensall,4)),1); % EMOS_tps(10,:,:)=nanmean(reshape(fr12km_EMOS_ensall,... % length(IMD_lat_points)*length(IMD_lon_points),... % size(fr12km_EMOS_ensall,3),size(fr12km_EMOS_ensall,4)),1); %% find 90th percentile including zeros. %to do percentile of all points with rain, just set zeros to NaNs IMD_tps_nozeros=IMD_tps; IMD_tps_nozeros(IMD_tps_nozeros==0)=NaN; IMD_90_rain=prctile(IMD_tps_nozeros,80,2); %repeat for IMERG %to do percentile of all points with rain, just set zeros to NaNs IMERG_tps_nozeros=IMERG_tps; IMERG_tps_nozeros(IMERG_tps_nozeros==0)=NaN; IMERG_90_rain=prctile(IMERG_tps_nozeros,75,2); %Answer the following question. At each timestep, how many ensemble members %are above this threshold? Fcst_thxced_rain=zeros(9,334); EMOS_thxced_rain=zeros(9,334); UQM_thxced_rain=zeros(9,334); for loc_c=1:9 Fcst_rain_exceed_mat=NaN(23,334); EMOS_rain_exceed_mat=NaN(23,334); UQM_rain_exceed_mat=NaN(23,334); Fcst_rain_exceed_mat(squeeze(Fcst_tps(loc_c,:,:))>IMERG_90_rain(loc_c))=1; Fcst_thxced_rain(loc_c,:)=nansum(Fcst_rain_exceed_mat,1)./23; EMOS_rain_exceed_mat(squeeze(EMOS_tps(loc_c,:,:))>IMERG_90_rain(loc_c))=1; EMOS_thxced_rain(loc_c,:)=nansum(EMOS_rain_exceed_mat,1)./23; UQM_rain_exceed_mat(squeeze(UQM_tps(loc_c,:,:))>IMERG_90_rain(loc_c))=1; UQM_thxced_rain(loc_c,:)=nansum(UQM_rain_exceed_mat,1)./23; end % Create matrix of same dimensions with 1 or 0's for if rainfall exceeded %threshold or not. IMERG_thxceed_mat=zeros(size(IMERG_tps)); for loc_IMD=1:9 A=IMD_tps(loc_IMD,:); IMERG_thxceed_mat(loc_IMD,:)=A>IMERG_90_rain(loc_IMD); end % % % loop to find how often the forecast is correct for each bin. % hist_edges=[0.0001 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1]; % Fcst_rd_mat=zeros(9,10); % UQM_rd_mat=zeros(9,10); % EMOS_rd_mat=zeros(9,10); % for xi=1:9 % fcst_X=Fcst_thxced_rain(xi,:); % UQM_X=UQM_thxced_rain(xi,:); % EMOS_X=EMOS_thxced_rain(xi,:); % for yj=1:10 % Fcst_rd_mat(xi,yj)=sum(IMD_thxceed_mat(xi,fcst_X>=hist_edges(yj) & fcst_X=hist_edges(yj) & fcst_X=hist_edges(yj) & UQM_X=hist_edges(yj) & UQM_X=hist_edges(yj) & EMOS_X=hist_edges(yj) & EMOS_X