latStart_string=strtrim(string(lat_start,format='(f5.1)'),2)
latEnd_string=strtrim(string(lat_end,format='(f5.1)'),2)
lonStart_string=strtrim(string(lon_start,format='(f6.1)'),2)
lonEnd_string=strtrim(string(lon_end,format='(f6.1)'),2)

dn_var=[1]
dn_string='Day'

dbz_thresholds=[-25.,-10.,0.,10.,30.]

;for zzz=0,n_elements(location_lon_vector)-1 do begin

;location_lat=location_lat_vector[zzz] & location_lon=location_lon_vector[zzz]
; ------------------------------------------------------------- ;
;location_lat = 0.5 * (lat_start + lat_end)
;location_lon = 0.5 * (lon_start + lon_end)
;dlat = lat_end - location_lat
;dlon = lon_end - location_lon
; ------------------------------------------------------------- ;

location_lat = 0.5 * (lat_start + lat_end)
dlat = lat_end - location_lat

if lon_start ge 0 then begin
   colon_start = lon_start
endif else begin
   colon_start = lon_start + 360.0
endelse

if lon_end ge 0 then begin
   colon_end = lon_end
endif else begin
   colon_end = lon_end + 360.0
endelse

colocation_lon = 0.5 * (colon_start + colon_end)

location_lon=colocation_lon
if colocation_lon gt 180. then location_lon=colocation_lon-360.

dlon = colon_end - colocation_lon

; ---------- ;

colon_vector=lon_vector
result=where(colon_vector lt 0,count)
if count gt 0 then colon_vector[result]=360.0+colon_vector[result]

; --------------------------------------------------------------- ;

i=0 
while lat_vector[i] le location_lat-dlat and i lt n_elements(lat_vector)-2 do i=i+1 & lat_start_index=i
while lat_vector[i] lt location_lat+dlat and i lt n_elements(lat_vector)-2 do i=i+1 & lat_end_index=i


i=0
while lon_vector[i] le location_lon-dlon and i lt n_elements(lon_vector)-2 do i=i+1 & lon_start_index=i
while lon_vector[i] lt location_lon+dlon and i lt n_elements(lon_vector)-2 do i=i+1 & lon_end_index=i

;i=0
;while colon_vector[i] le colocation_lon-dlon and i lt n_elements(colon_vector)-2 do i=i+1 & lon_start_index=i
;while colon_vector[i] lt colocation_lon+dlon and i lt n_elements(colon_vector)-2 do i=i+1 & lon_end_index=i

; --------------------------------------------------------------- ;
;if location_lat eq 0. and location_lon eq 0. then begin
;    lat_start_index=0. & lat_end_index=n_elements(lat_vector)-1
;    lon_start_index=0. & lon_end_index=n_elements(lon_vector)-1
;endif
; --------------------------------------------------------------- ;

vertical_frequency=fltarr(n_elements(height_vector),n_elements(dbz_thresholds))
;Geoprof_Occurrence(Cloudsat_Layer_Top_Type, MODIS_scene_variability, MODIS_scene_characterization, daytime_flag, longitude, latitude)

; a 2-d array that stores the stats of the types

cpr_modis_type_histo=fltarr(6,10)

for j=0,n_elements(height_vector)-1 do begin
    ;total_CPR_Occ_Grid = total(float(CPR_Occ_Grid[lat_start_index:lat_end_index, lon_start_index:lon_end_index, j, *, dn_var]))
    total_CPR_Occ_Grid = total(float(CPR_Occ_Grid[lat_start_index:lat_end_index, lon_start_index:lon_end_index, j, *, *]))
    ;print, "j, total_CPR_Occ_Grid: ", j, total_CPR_Occ_Grid
for k=0,n_elements(dbz_thresholds)-1 do begin

    i=0
    while dbz_vector[i] lt dbz_thresholds[k] and i lt n_elements(dbz_vector)-2 do i=i+1 & dbz_end_index=i
    
    if total_CPR_Occ_Grid ne 0 then begin
        ;total_CPR_Occ_Grid_K = total(float(CPR_Occ_Grid[lat_start_index:lat_end_index, lon_start_index:lon_end_index, j, 1:dbz_end_index, dn_var]))
        total_CPR_Occ_Grid_K = total(float(CPR_Occ_Grid[lat_start_index:lat_end_index, lon_start_index:lon_end_index, j, 1:dbz_end_index, *]))
        vertical_frequency[j,k]=total_CPR_Occ_Grid_K / total_CPR_Occ_Grid
    endif else begin    ; if total(float(CPR_Occ_Grid[lat_start_indx:lat_end_indx,jj,k,*,*])) ne 0 then begin
       vertical_frequency[j,k]=0.0
    endelse ; if total(float(CPR_Occ_Grid[lat_start_indx:lat_end_indx,jj,k,*,*])) ne 0 then begin

endfor
endfor       ; for j=0,n_elements(height_vector)-1

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; geoprof occurrence
;
; limit the scene var to highly uniform
for jj=0,5 do begin ; cpr echo tope
for kk=0,9 do begin ; modis scene type

    if total(float(geoprof_occ_grid[lat_start_index:lat_end_index,lon_start_index:lon_end_index,*,kk,1,jj])) ne 0 then begin
        cpr_modis_type_histo[jj,kk]=total(float(geoprof_occ_grid[lat_start_index:lat_end_index,lon_start_index:lon_end_index,*,kk,1,jj]));/$
        ;                  total(float(geoprof_occ_grid[lat_start_index:lat_end_index,lon_start_index:lon_end_index,*,*,1,jj]))
    endif else begin    ; 
        cpr_modis_type_histo[jj,kk]=0.0
    endelse ; 


endfor  ; for kk=0,10 do begin
endfor  ; for jj=0,5 do begin