; $Id: dm_plot_trajectory.pro,v 1.23 2008/03/19 20:46:50 ymqiu Exp $
;#######################################################################
;
; NAME:
;  dm_plot_trajectory
;
; PURPOSE:
;  plot trajectory
;
; CATEGORY:
;  dcs_mslice
;
; AUTHOR:
;  Yiming Qiu
;  NIST Center for Neutron Research
;  100 Bureau Drive, Gaithersburg, MD 20899-6102
;  United States
;  yiming.qiu@nist.gov
;  March, 2008
;
; LICENSE:
;  The software in this file is written by an employee of
;  National Institute of Standards and Technology
;  as part of the DAVE software project.
;
;  The DAVE software package is not subject to copyright protection
;  and is in the public domain. It should be considered as an
;  experimental neutron scattering data reduction, visualization, and
;  analysis system. As such, the authors assume no responsibility
;  whatsoever for its use, and make no guarantees, expressed or
;  implied, about its quality, reliability, or any other
;  characteristic. The use of certain trade names or commercial
;  products does not imply any endorsement of a particular product,
;  nor does it imply that the named product is necessarily the best
;  product for the stated purpose. We would appreciate acknowledgment
;  if the DAVE software is used or if the code in this file is
;  included in another product.
;
;#######################################################################

pro dcs_mslice::dm_plot_trajectory,overplot=overplot,debug=debug
    if keyword_set(debug) then current=systime(/sec)
    if (not ptr_valid(self.projStrPtr) ) or self.error then begin
       ok = dialog_message('Please calculate projection first.',/error,dialog_parent=self.tlb)
       return
    endif

    if self.ftype eq 0 then none_id = [5,7,5,6] else none_id = [4,7,5,6]

    if (self.traj_x[self.samp_typ] eq none_id[self.samp_typ]) or (self.traj_y[self.samp_typ] eq none_id[self.samp_typ]) then begin
       ok = dialog_message('Invalid axis.',/error,dialog_parent=self.tlb)
       return
    endif

    ;select data
    self->dm_choose_trajdata,self.traj_x[self.samp_typ],data=xdat,label=xtit
    self->dm_choose_trajdata,self.traj_y[self.samp_typ],data=ydat,label=ytit
    if xtit eq 'E (meV)' or xtit eq 'Energy' then begin
       if self.samp_typ eq 3 then uniq_E=xdat[*,0,0] $
       else uniq_E=xdat[*,0]
    endif
    if ytit eq 'E (meV)' or ytit eq 'Energy' then begin
       if self.samp_typ eq 3 then uniq_E=ydat[*,0,0] $
       else uniq_E=ydat[*,0]
    endif

    if (self.traj_z[self.samp_typ] ne none_id[self.samp_typ]) and (total(finite(self.traj_zran)) ne 0) then begin
       self->dm_choose_trajdata,self.traj_z[self.samp_typ],range=self.traj_zran,index=index,empty=empty
       if n_elements(xdat) ne n_elements(index) then begin
          xdat = xdat[index]
          ydat = ydat[index]
       endif
       if empty then begin
          ok = dialog_message('The specifiled z range covers no data. The entire possible range is used.',/error,dialog_parent=self.tlb)
       endif
       index = -1
    endif

    ;xrange and yrange
    n_xran = total(finite(self.traj_xran))
    n_yran = total(finite(self.traj_yran))
    if n_xran eq 2 then begin
       if self.traj_xran[1] gt self.traj_xran[0] then begin
          index = where((xdat le self.traj_xran[1]) and (xdat ge self.traj_xran[0]),count)
          if count gt 0 then begin
             xdat  = xdat[index]
             ydat  = ydat[index]
             index = -1
          endif
       endif
    endif
    if n_yran eq 2 then begin
       if self.traj_yran[1] gt self.traj_yran[0] then begin
          index = where((ydat le self.traj_yran[1]) and (ydat ge self.traj_yran[0]),count)
          if count gt 0 then begin
             xdat  = xdat[index]
             ydat  = ydat[index]
             index = -1
          endif
       endif
    endif
    if n_xran lt 2 then begin
       if finite(self.traj_xran[0]) then begin
          index = where(xdat ge self.traj_xran[0],count)
          if count gt 0 then begin
             xdat  = xdat[index]
             ydat  = ydat[index]
             index = -1
          endif
       endif
       if finite(self.traj_xran[1]) then begin
          index = where(xdat le self.traj_xran[1],count)
          if count gt 0 then begin
             xdat  = xdat[index]
             ydat  = ydat[index]
             index = -1
          endif
       endif
    endif
    if n_yran lt 2 then begin
       if finite(self.traj_yran[0]) then begin
          index = where(ydat ge self.traj_yran[0],count)
          if count gt 0 then begin
             xdat  = xdat[index]
             ydat  = ydat[index]
             index = -1
          endif
       endif
       if finite(self.traj_yran[1]) then begin
          index = where(ydat le self.traj_yran[1],count)
          if count gt 0 then begin
             xdat  = xdat[index]
             ydat  = ydat[index]
             index = -1
          endif
       endif
    endif

    if n_elements(xdat) gt 1e5 then begin
       xmax = max(xdat,min=xmin)
       ymax = max(ydat,min=ymin)
       nogrid = 0
       if (xtit eq 'E (meV)') or (xtit eq 'Energy') then begin
          nogrid = 1
          index  = where((uniq_E ge xmin) and (uniq_E le xmax),count)
          if count ne 0 then uniq_E=uniq_E[index]
          index  = -1
       endif
       if (ytit eq 'E (meV)') or (ytit eq 'Energy') then begin
          nogrid = 2
          index  = where((uniq_E ge ymin) and (uniq_E le ymax),count)
          if count ne 0 then uniq_E=uniq_E[index]
          index  = -1
       endif
       xstep = (xmax-xmin)/300.
       ystep = (ymax-ymin)/300.
       if xstep eq 0 then xstep=0.1
       if ystep eq 0 then ystep=0.1
       zdat  = bytarr(n_elements(xdat))+1b
       if nogrid eq 1 then extra={uniq_xval:uniq_E} $
       else if nogrid eq 2 then extra={uniq_yval:uniq_E}
       dm_stepgrid_bin,xstep,ystep,xdat,ydat,zdat,debug=debug,_extra=extra
       index = where(finite(temporary(zdat)),count)
       if count ne 0 then begin
          n_x   = n_elements(xdat)
          xdat  = xdat[index mod n_x]
          ydat  = ydat[fix(index/n_x,type=3)]
          index = -1
       endif
    endif
    if keyword_set(debug) then help,xdat

    ;extra info in the title
    if self.ftype eq 0 then $
       title = dm_to_string(fix(sqrt(81.8042/self.eief)*100.+0.5)/100.)+string('c5'XB) $
    else $
       title = (['Ei=','Ef='])[self.instrgeom]+dm_to_string(fix(self.eief*1000.+0.5,type=3)/1000.)+' meV'
    switch self.samp_typ of
        1:  begin
            title = title+', psi='+dm_to_string(self.orie_psi)+string('b0'XB)
            break
            end
        2:
        3:  begin
            title = title+', psi=['+dm_to_string(self.a2[0])+string('b0'XB)+$
                    ':'+dm_to_string(self.a2[1])+string('b0'XB)+':'+$
                    dm_to_string(self.a2[2])+string('b0'XB)+']'
            break
            end
    endswitch
    tmp = ''
    if (self.traj_z[self.samp_typ] ne none_id[self.samp_typ]) and (total(finite(self.traj_zran)) eq 2) then begin
       case self.samp_typ of
        0:  types = self.view_type
        1:  types = self.traj_type
        2:  types = [self.traj_type[0:4],self.traj_type[7]]
        3:  types = [self.traj_type[0:5],self.traj_type[7]]
        else:
       endcase
       tmp = types[self.traj_z[self.samp_typ]]
       if (strlowcase(tmp) eq 'e') or (strlowcase(tmp) eq 'energy') or (strlowcase(tmp) eq 'e (mev)') $
          then tmp='E'
       title = title+', '+tmp+'=['+dm_to_string(self.traj_zran[0])+','+dm_to_string(self.traj_zran[1])+']'
       if tmp eq 'E' then title=title+'meV'
    endif
    if (self.samp_typ ge 2) and (tmp ne 'E') then begin
       case self.e_bin[self.samp_typ] of
          0:  title = title
          1:  title = title+', Elastic'
          2:  title = title+', E=['+dm_to_string(self.e_range[0])+','+dm_to_string(self.e_range[1])+']meV'
       endcase
    endif

    ;plot the data
    isnew = 0b
    if self.useitool then begin
       ok = dialog_message('Itools plotting is not implemented yet. Soooorry!',dialog_parent=self.tlb)
       return
    endif else begin
       obj_name = 'dm_plot'
       if not obj_valid(self.current) then begin
          self.current=obj_new(obj_name,xdat,ydat,xtit=xtit,ytit=ytit,/no_copy,/hidelegend,/showxgrid,$
                /showygrid,parentobj=self,linestyle='no line',psym='dot',color='green',xran=self.traj_xran,$
                yran=self.traj_yran,layer=5,title=title,legend='trajectory')
          isnew = 1b & self.traj_layer=0
       endif else begin
          self.current->getproperty,keep=keep
          if keep eq 1 then begin
             isnew = 1b
             self.current=obj_new(obj_name,linestyle='no line',psym='dot',parentobj=self,/hidelegend)
             self.traj_layer = 0
          endif
          if keyword_set(overplot) then begin
             if not isnew then self.traj_layer=self.traj_layer+1
             self.current->add_plot,xdat,ydat+0.005*self.traj_layer*(-1)^(self.traj_layer),/no_copy,$
                psym='dot',linestyle='no line',color=(['green','red','magenta','cyan','blue'])[(self.traj_layer mod 5)],$
                layer=5+self.traj_layer,legend='trajectory'
          endif else $
             self.current->setproperty,xdat=xdat,ydat=ydat,xtit=xtit,ytit=ytit,/no_copy,ylog=0,$
                xran=self.traj_xran,yran=self.traj_yran,color='green',psym='dot',linestyle='no line',$
                /hidelegend,title=title,/nodraw,layer=5,/showxgrid,/showygrid,legend='trajectory'
       endelse
       ;add powder lines in diffuse scattering mode
       if ptr_valid(self.powd_line) and (isnew or (not keyword_set(overplot))) then begin
          tmp = size(*self.powd_line)
          tmp = tmp[tmp[0]]
          if self.traj_x[self.samp_typ] eq 3 and self.traj_y[self.samp_typ] eq 4 then begin
             ix = 0 & iy = 1
          endif else if self.traj_x[self.samp_typ] eq 4 and self.traj_y[self.samp_typ] eq 3 then begin
             ix = 1 & iy = 0
          endif
          if n_elements(ix) ne 0 then $
             for i=0,tmp-1 do $
                 self.current->add_plot,(*self.powd_line)[*,ix,i],(*self.powd_line)[*,iy,i],$
                        psym='no symbol',color='yellow',linestyle='solid',layer=i,$
                        legend=(['Al','Cu','Stainless Steel'])[(self.powd_typ-1)>0]+' powder line'
       endif
       self.current->draw
    endelse
    if keyword_set(overplot) and (not isnew)then $
       self->dm_add_plotwindow,type=2,title=strtrim((strsplit(xtit,'()',/extract))[0],2)+' vs. '+$
           strtrim((strsplit(ytit,'()',/extract))[0],2)+'[overplot]'  $
    else $
       self->dm_add_plotwindow,type=2,title=strtrim((strsplit(xtit,'()',/extract))[0],2)+' vs. '+$
           strtrim((strsplit(ytit,'()',/extract))[0],2)
    if keyword_set(debug) then print,'plot_trajectory finished in ',systime(/sec)-current,' sec.'
end