pro tile_check_1c,user_name=user_name,trial=trial,orbit_num=orbit_num,fix_sigma=fix_sigma,fix_c=fix_c,n_clicks=n_clicks,x=x,y=y
  common curvefit_trace,rad,input_c,input_sigma,input_a,input_r,coeff_record
  coeff_record=0
  input_c=39
  input_sigma=0.608d
  input_r=30
  input_a=25

  ;Clear the window, and get the rainbow+white color bar
  window,0,xpos=100,ypos=0,xsize=1800,ysize=1200
  erase
  device,decomposed=0
  loadct,39

  window,1,xpos=1921,ypos=0,xsize=1280,ysize=1024
  erase
  device,decomposed=0
  loadct,39

  window,2,xpos=2121,ypos=0,xsize=1100,ysize=1100
  erase
  device,decomposed=0
  loadct,39

  ;zoom is number of display pixels per tile
  zoom=5

  ;these ofs's let you shift the image on the screen
  xofs=-200
  yofs=0

  ;Load the data
  if n_elements(trial) eq 0 then trial=85
  if n_elements(user_name) eq 0 then user_name="jeppesen"
  if n_elements(orbit_num) eq 0 then orbit_num=0
;  restore,get_level_1c_filename(trial_num=trial,orbit_num,user_name=user_name)
;  print,get_level_1c_filename(trial_num=trial,orbit_num,user_name=user_name)
;  restore,get_debug_filename("tile",trial_num=trial,orbit_num,user_name=user_name)
;  print,get_debug_filename("tile",trial_num=trial,orbit_num,user_name=user_name)
;  restore,get_debug_filename("cloud_initial_guess",trial_num=trial,orbit_num,user_name=user_name)
;  print,get_debug_filename("cloud_initial_guess",trial_num=trial,orbit_num,user_name=user_name)

  ;Plot the image
  wset,0
  s=size(*(level_1c_retrievals.scale_height_ratio),/dimensions)
  tvscl,/nan,rebin(*(level_1c_retrievals.albedo),s[0]*zoom,s[1]*zoom,/sample),xofs*zoom,yofs*zoom
  if n_elements(n_clicks) eq 0 then n_clicks=1

  ;Get 10 mouse clicks
  for i=1,n_clicks do begin
    wset,0
    junk=temporary(coeff_record)
    ;Get the mouse click
    if i eq 1 and n_elements(x) eq 0 then begin
      cursor,x,y,/device,/up
      x=x/zoom-xofs
      y=y/zoom-yofs
    end
    print,x,y

    ;X marks the spot
    plots,([x-5,x+5]+xofs)*zoom+zoom/2,([y-5,y+5]+yofs)*zoom+zoom/2,/device,color=255
    plots,([x-5,x+5]+xofs)*zoom+zoom/2,([y+5,y-5]+yofs)*zoom+zoom/2,/device,color=255

    wset,1

    print,string(x,y,'(%"this_tiles.ssa[%03d,%03d,*]:")')
    print,string(tiles.ssa[x,y,*],format='(%"%20.6f")')
    print,string(x,y,'(%"this_tiles.rad[%03d,%03d,*] (G):")')
    print,string(tiles.rad[x,y,*],format='(%"%20.6f")')

    print,string(x,y,(*(level_1c_retrievals.cloud_no_phase))[x,y],'(%"actual_albedo[%03d,%03d]: %20.6fG")')
    input_a=(*(level_1c_retrievals.cloud_no_phase))[x,y]
    print,string(x,y,(*(level_1c_retrievals.albedo))[x,y],'(%"albedo[%03d,%03d]: %20.6fG")')
    print,string(x,y,(*(level_1c_retrievals.particle_radius))[x,y],'(%"radius[%03d,%03d]: %20.6fnm")')
    print,string(x,y,(*(level_1c_retrievals.ozone_col_density))[x,y],'(%"ozone[%03d,%03d]: %20.6fPeta")')
    print,string(x,y,(*(level_1c_retrievals.scale_height_ratio))[x,y],'(%"hratio[%03d,%03d]: %20.6f")')
    cm=(*(level_1c_retrievals.correlation))[x,y,*]
    c=(*(level_1c_retrievals.ozone_col_density_unc))[x,y,*]
    h=(*(level_1c_retrievals.scale_height_ratio_unc))[x,y,*]
    a=(*(level_1c_retrievals.albedo_unc))[x,y,*]
    r=(*(level_1c_retrievals.particle_radius_unc))[x,y,*]
    corr=fltarr(4,4)
    corr[[0,5,10,15]]=1
    corr[0,1:3]=cm[0:2]
    corr[1,2:3]=cm[3:4]
    corr[2,3]=cm[5]
    corr[1:3,0]=cm[0:2]
    corr[2:3,1]=cm[3:4]
    corr[3,2]=cm[5]
    print,"stored correlation matrix"
    print,corr
    sig=[c,h,a,r]
    print,"stored sigmas"
    print,sig
    sigmtx=dblarr(4,4)
    for i_sigmtx=0,3 do for j_sigmtx=0,3 do begin
      sigmtx[i_sigmtx,j_sigmtx]=sig[i_sigmtx]*sig[j_sigmtx]
    end
    covar=corr*sigmtx
    print,"stored covariance matrix"
    print,covar

    print,string(x,y,initial_guess[x,y,2],'(%"ig_albedo[%03d,%03d]: %20.6fG")')
    print,string(x,y,initial_guess[x,y,3],'(%"radius[%03d,%03d]: %20.6fnm")')
    print,string(x,y,initial_guess[x,y,0],'(%"ozone[%03d,%03d]: %20.6fPeta")')
    print,string(x,y,initial_guess[x,y,1],'(%"hratio[%03d,%03d]: %20.6f")')

    ;Draw a plot of input SSA (xaxis) vs input radiance (yaxis) for all layers in this tile
    plot,tiles.ssa[x,y,*],tiles.rad[x,y,*],charsize=2,xrange=[0,180],yrange=[0,500],/noerase

    ;plot,(*(level_1c_retrievals.scale_height_ratio))[*,y],/noerase,charsize=2

    ;Trace the curve fit
    w=where(finite(tiles.rad[x,y,*]),nw)
    if nw gt 0 then rad=tiles.rad[x,y,w] else rad=tiles.rad[x,y,*]
    t={ BINX:tiles.binx, $
   BINY      :tiles.BINY, $
   STARTPIXX :tiles.STARTPIXX, $
   STARTPIXY :tiles.STARTPIXY , $
   PIXX      :tiles.PIXX      , $
   PIXY      :tiles.PIXY      , $
   N_X       :tiles.N_X       , $
   N_XT      :tiles.N_XT      , $
   N_Y       :tiles.N_Y       , $
   N_YT      :tiles.N_YT      , $
   N_GOOD    :tiles.N_GOOD[x,y,*], $
   SSA       :tiles.SSA[x,y,*]      , $
   SZA       :tiles.SZA[x,y,*]       , $
   SVA       :tiles.SVA[x,y,*]       , $
   CNP       :tiles.CNP[x,y,*]       , $
   RAD       :tiles.RAD[x,y,*]       , $
   RUC       :tiles.RUC[x,y,*]        $
   }
    c=findgen(10)
      ;Stop tracing the curve fit
    junk=temporary(rad)
    for j=0,n_elements_c-1 do begin
;      result=calc_cloud_parm(t,1,firstc=initial_guess[x,y,0],/verbose,covar=covar,fix_c=fix_c,fix_sigma=fix_sigma)
    end
;    model_difference=(t.rad[w]-calc_obsbrt_cloud(t.ssa[w], [input_c,input_sigma,input_a,input_r], sva=t.sva[w],sza=t.sza[w]))
    model_norm_difference=model_difference/t.ruc[w]
    model_variance=model_norm_difference^2
    model_chisq=total(model_variance)
    print,"correct chisq: ",model_chisq
    print,gammq((nw-4)/2,model_chisq/2)
;    difference=(t.rad[w]-calc_obsbrt_cloud(t.ssa[w], [result.ret_o3cden,result.ret_hratio,result.ret_clda,result.ret_cldr], sva=t.sva[w],sza=t.sza[w]))
    norm_difference=difference/t.ruc[w]
    variance=norm_difference^2
    chisq=total(variance)
    print,"fit chisq: ",chisq
    print,gammq((nw-4)/2,chisq/2)
    print,"covariance:"
    covar=fltarr(4,4)
    covar[[0,5,10,15]]=1
    covar[0,1:3]=result.corr[0,0,0:2]
    covar[1,2:3]=result.corr[0,0,3:4]
    covar[2,3]=result.corr[0,0,5]
    covar[1:3,0]=result.corr[0,0,0:2]
    covar[2:3,1]=result.corr[0,0,3:4]
    covar[3,2]=result.corr[0,0,5]
    sig=[result.sig_o3cden,result.sig_hratio,result.sig_clda,result.sig_cldr]
    sigmtx=dblarr(4,4)
    for i_sigmtx=0,3 do for j_sigmtx=0,3 do begin
      sigmtx[i_sigmtx,j_sigmtx]=sig[i_sigmtx]*sig[j_sigmtx]
    end
    covar=covar*sigmtx

    print,covar
    print,"sigma: ",sig
    sigmtx=dblarr(4,4)
    for i_sigmtx=0,3 do for j_sigmtx=0,3 do begin
      sigmtx[i_sigmtx,j_sigmtx]=sig[i_sigmtx]*sig[j_sigmtx]
    end
    corr=covar/sigmtx
    print,"correlation:"
    print,corr
    wset,2
;    plot,coeff_record[0,*],coeff_record[2,*],xtitle="Column ozone density (P)",ytitle="Cloud albedo (G)",xrange=[0,50],yrange=[-10,50],psym=2,/noerase,color=i*(254-60)+60
;    print,coeff_record
;    parameter space
     w=where(finite(tiles.rad[x,y,*]),nw)
     ssa=tiles.ssa[x,y,w]
     sva=tiles.sva[x,y,w]
     sza=tiles.sza[x,y,w]
     rad_notrace=tiles.rad[x,y,w]
    if n_elements(a) gt 0 then junk=temporary(a)
    junk=0
    for i_pspace=0,9 do begin
      print,i_pspace
      parameter_space_step,a,rad_notrace,ssa,sva,sza
      print,n_elements(a)
      print,(size(a,/dimensions))[0]
      save,a,filename="parameter_space.sav"
      if i_pspace ne 10 then tvscl,congrid(alog(a),1025,1025) else tvscl,alog(a)
      plots,/device,coeff_record[2,*]*1025d/40d,coeff_record[3,*]*1025d/100d,psym=2,color=255
    end
;    tvscl,alog(a)
  end
end