; $Id: gaussian98.pro,v 1.8 2005/07/21 14:08:45 craig_nist Exp $
; 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 of if the code in this file is
;  included in another product.
;


;   Read Gaussian files


;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gausscleanup,tlb
  widget_control,tlb,get_uvalue = pState

if ptr_valid(((*pState).atomorig).massptr) then ptr_free,((*pState).atomorig).massptr
if ptr_valid(((*pState).atomnew).newmassptr) then ptr_free,((*pState).atomnew).newmassptr
if ptr_valid(((*pState).atomorig).listptr) then ptr_free,((*pState).atomorig).listptr
if ptr_valid(((*pState).atomnew).newlistptr) then ptr_free,((*pState).atomnew).newlistptr
if ptr_valid(((*pState).atomorig).numberptr) then ptr_free,((*pState).atomorig).numberptr
if ptr_valid(((*pState).atomnew).newnumberptr) then ptr_free,((*pState).atomnew).newnumberptr
if ptr_valid(((*pState).atomorig).ixsectionptr) then ptr_free,((*pState).atomorig).ixsectionptr
if ptr_valid(((*pState).atomnew).newixsectionptr) then ptr_free,((*pState).atomnew).newixsectionptr


                    if ptr_valid((*(*pState).readdataPtr).DatePtr) then ptr_free,(*(*pState).readdataPtr).DatePtr
                    if ptr_valid((*(*pState).readdataPtr).ScanPtr) then ptr_free,(*(*pState).readdataPtr).ScanPtr
                    if ptr_valid((*(*pState).readdataPtr).BasePtr) then ptr_free,(*(*pState).readdataPtr).BasePtr
                    if ptr_valid((*(*pState).readdataPtr).TypePtr) then ptr_free,(*(*pState).readdataPtr).TypePtr
                    if ptr_valid((*(*pState).readdataPtr).Collimationptr) then ptr_free,(*(*pState).readdataPtr).Collimationptr
                    if ptr_valid((*(*pState).readdataPtr).MosaicPtr) then ptr_free,(*(*pState).readdataPtr).MosaicPtr
                    if ptr_valid((*(*pState).readdataPtr).XtalOrientPtr) then ptr_free,(*(*pState).readdataPtr).XtalOrientPtr
                    if ptr_valid((*(*pState).readdataPtr).LatticeAPtr) then ptr_free,(*(*pState).readdataPtr).LatticeAPtr
                    if ptr_valid((*(*pState).readdataPtr).LatticeBPtr) then ptr_free,(*(*pState).readdataPtr).LatticeBPtr
                    if ptr_valid((*(*pState).readdataPtr).LatticeCPtr) then ptr_free,(*(*pState).readdataPtr).LatticeCPtr
                    if ptr_valid((*(*pState).readdataPtr).LatticeAlPtr) then ptr_free,(*(*pState).readdataPtr).LatticeAlPtr
                   if ptr_valid( (*(*pState).readdataPtr).LatticeBePtr) then ptr_free,(*(*pState).readdataPtr).LatticeBePtr
                    if ptr_valid((*(*pState).readdataPtr).LatticeGaPtr) then ptr_free,(*(*pState).readdataPtr).LatticeGaPtr
                    if ptr_valid((*(*pState).readdataPtr).E_ZEROPtr) then ptr_free,(*(*pState).readdataPtr).E_ZEROPtr
                   if ptr_valid( (*(*pState).readdataPtr).DELTA_EPtr) then ptr_free,(*(*pState).readdataPtr).DELTA_EPtr
                    if ptr_valid((*(*pState).readdataPtr).E_FINALPtr) then ptr_free,(*(*pState).readdataPtr).E_FINALPtr
                    if ptr_valid((*(*pState).readdataPtr).MON_D_SPACEPtr) then ptr_free,(*(*pState).readdataPtr).MON_D_SPACEPtr
                    if ptr_valid((*(*pState).readdataPtr).AN_D_SPACEPtr) then ptr_free,(*(*pState).readdataPtr).AN_D_SPACEPtr
                    if ptr_valid((*(*pState).readdataPtr).T_startPtr) then ptr_free,(*(*pState).readdataPtr).T_startPtr
                    if ptr_valid((*(*pState).readdataPtr).T_INCPtr) then ptr_free,(*(*pState).readdataPtr).T_INCPtr
                    if ptr_valid((*(*pState).readdataPtr).HH_Ptr) then ptr_free,(*(*pState).readdataPtr).HH_Ptr
                    if ptr_valid((*(*pState).readdataPtr).KK_Ptr) then ptr_free,(*(*pState).readdataPtr).KK_Ptr
                    if ptr_valid((*(*pState).readdataPtr).LL_Ptr) then ptr_free,(*(*pState).readdataPtr).LL_Ptr
                    if ptr_valid((*(*pState).readdataPtr).dH_Ptr) then ptr_free,(*(*pState).readdataPtr).dH_Ptr
                    if ptr_valid((*(*pState).readdataPtr).dK_Ptr) then ptr_free,(*(*pState).readdataPtr).dK_Ptr
                    if ptr_valid((*(*pState).readdataPtr).dL_Ptr) then ptr_free,(*(*pState).readdataPtr).dL_Ptr
                    if ptr_valid((*(*pState).readdataPtr).HFIELD_Ptr) then ptr_free,(*(*pState).readdataPtr).HFIELD_Ptr
                    if ptr_valid((*(*pState).readdataPtr).MonitorPtr) then ptr_free,(*(*pState).readdataPtr).MonitorPtr
                    if ptr_valid((*(*pState).readdataPtr).PrefactorPtr) then ptr_free,(*(*pState).readdataPtr).PrefactorPtr
                    if ptr_valid((*(*pState).readdataPtr).QxPtr) then ptr_free,(*(*pState).readdataPtr).QxPtr
                    if ptr_valid((*(*pState).readdataPtr).QyPtr) then ptr_free,(*(*pState).readdataPtr).QyPtr
                    if ptr_valid((*(*pState).readdataPtr).QzPtr) then ptr_free,(*(*pState).readdataPtr).QzPtr
                    if ptr_valid((*(*pState).readdataPtr).EnergyPtr) then ptr_free,(*(*pState).readdataPtr).EnergyPtr
                    if ptr_valid((*(*pState).readdataPtr).TempPtr) then ptr_free,(*(*pState).readdataPtr).TempPtr
                    if ptr_valid((*(*pState).readdataPtr).TimePtr) then ptr_free,(*(*pState).readdataPtr).TimePtr
                    if ptr_valid((*(*pState).readdataPtr).TotalIntPtr) then ptr_free,(*(*pState).readdataPtr).TotalIntPtr
           if ptr_valid(((*pState).readdataPtr)) then ptr_free,((*pState).readdataPtr)

  if ptr_valid((*pstate).initialorientationptr) then ptr_free,(*pstate).initialorientationptr
  if ptr_valid((*pstate).isotopeinfoptr) then ptr_free,(*pstate).isotopeinfoptr
  if ptr_valid((*pstate).eigenyptr) then ptr_free,(*pstate).eigenyptr
  if ptr_valid((*pstate).spectrumptr) then ptr_free,(*pstate).spectrumptr
;  if ptr_valid((*pstate).inptr) then ptr_free,(*pstate).inptr
  if ptr_valid((*pState).symmodeptr)then ptr_free,(*pState).symmodeptr
  if ptr_valid((*pState).freqptr)then ptr_free,(*pState).freqptr
  if ptr_valid((*pState).origfreqptr)then ptr_free,(*pState).origfreqptr
  if ptr_valid((*pState).redmassptr)then ptr_free,(*pState).redmassptr
  if ptr_valid((*pState).frccptr)then ptr_free,(*pState).frccptr
  if ptr_valid((*pState).iriptr)then ptr_free,(*pState).iriptr
  if ptr_valid((*pState).ramanptr)then ptr_free,(*pState).ramanptr
  if ptr_valid((*pState).depolarptr)then ptr_free,(*pState).depolarptr
  if ptr_valid((*pState).vectorptr)then ptr_free,(*pState).vectorptr

  if ptr_valid(pstate) then ptr_free,pstate

end
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gausshelp,event
;pdf_file = file_which('gaussian1.pdf',/include_current_dir)
pdf_file = !DAVE_AUXILIARY_DIR+'gaussian1.pdf'
void = launch_help(pdf_file,tlb = event.top)
return
end
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gaussabout,event
strout=["Turn your Gaussian98 or G03 frequency calculation into a useful comparison to FANS data!",$
       "brought to you by","Craig Brown (2002)","","craig.brown@nist.gov "]
void=dialog_message(dialog_parent=event.top,strout,/information)
return
end

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gaussquit,event
  widget_control,event.top,/destroy
end
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

pro gaussChooseFile,event
widget_control, event.top, get_uvalue = pState
      (*pState).gfile = DIALOG_PICKFILE(/READ,/MUST_EXIST,$
          FILTER='*.log')

if n_elements((*pState).gfile) gt 0 then begin
;widget_control,event.top,set_uvalue=pState,/no_copy
;clear all pointers to do with gaussian file to be '0'
        widget_control,event.top,get_uvalue=pState,/no_copy

  if ptr_valid(((*pState).atomorig).massptr) then *((*pState).atomorig).massptr=''
  if ptr_valid(((*pState).atomnew).newmassptr) then *((*pState).atomnew).newmassptr=''
  if ptr_valid(((*pState).atomorig).listptr) then *((*pState).atomorig).listptr=''
  if ptr_valid(((*pState).atomnew).newlistptr) then *((*pState).atomnew).newlistptr=''
  if ptr_valid(((*pState).atomorig).numberptr) then *((*pState).atomorig).numberptr=''
  if ptr_valid(((*pState).atomnew).newnumberptr) then *((*pState).atomnew).newnumberptr=''
  if ptr_valid(((*pState).atomorig).ixsectionptr) then *((*pState).atomorig).ixsectionptr=''
  if ptr_valid(((*pState).atomnew).newixsectionptr) then *((*pState).atomnew).newixsectionptr=''

    if ptr_valid((*pstate).eigenyptr) then *(*pstate).eigenyptr=''
    if ptr_valid((*pstate).spectrumptr) then *(*pstate).spectrumptr=''
;    if ptr_valid((*pstate).inptr) then *(*pstate).inptr=''
    if ptr_valid((*pState).symmodeptr)then *(*pState).symmodeptr=''
    if ptr_valid((*pState).freqptr)then *(*pState).freqptr=''
    if ptr_valid((*pState).origfreqptr)then *(*pState).origfreqptr=''
    if ptr_valid((*pState).redmassptr)then *(*pState).redmassptr=''
    if ptr_valid((*pState).frccptr)then *(*pState).frccptr=''
    if ptr_valid((*pState).iriptr)then *(*pState).iriptr=''
    if ptr_valid((*pState).ramanptr)then *(*pState).ramanptr=''
    if ptr_valid((*pState).depolarptr)then *(*pState).depolarptr=''
    if ptr_valid((*pState).vectorptr)then *(*pState).vectorptr=''
  (*pState).rotconstants=''
  (*pState).lowfrequencies=''
  (*pState).job=''                  ; job type array (strings)
  (*pState).dof=0                           ; degrees of freedom  (int)
  (*pState).maxforce=''              ; max force obtained, and threshold.
  (*pState).maxforceflag=''                 ; did it meet the criteria
  (*pState).rmsforce=''              ; RMS force obtained, and threshold.
  (*pState).rmsforceflag=''                 ; did it meet the criteria
  (*pState).maxdisp=''               ; max disp obtained, and threshold.
  (*pState).maxdispflag=''                  ; did it meet the criteria
  (*pState).rmsdisp=''               ; RMS disp obtained, and threshold.
  (*pState).rmsdispflag=''                  ; did it meet the criteria
  (*pState).numberatoms=0                   ; how many atoms in the calculation
widget_control,event.top,set_uvalue=pState,/no_copy
   gaussReadGaussian,event
endif
widget_control,event.top,set_uvalue=pState,/no_copy
end
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gaussdump,event
widget_control, event.top, get_uvalue = pState
dump=(*pstate).dump
if (dump eq 0) then (*pstate).dump=1
if (dump eq 1) then (*pstate).dump=0
widget_control, event.top, set_uvalue = pState
gaussplotcalc,event
end
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gaussdumpdave,event
widget_control, event.top, get_uvalue = pState
dumpdave=(*pstate).dumpdave
if (dumpdave eq 0) then (*pstate).dumpdave=1
if (dumpdave eq 1) then (*pstate).dumpdave=0
widget_control, event.top, set_uvalue = pState
gaussplotcalc,event
end
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gaussChooseData,event
  widget_control, event.top, get_uvalue = pState
    (*pState).dfile = DIALOG_PICKFILE(/READ,/MUST_EXIST,$
          FILTER='*.*')

if n_elements((*pState).dfile) gt 0 then begin
                      if ptr_valid((*(*pState).readdataPtr).DatePtr) then *(*(*pState).readdataPtr).DatePtr=''
                    if ptr_valid((*(*pState).readdataPtr).ScanPtr) then *(*(*pState).readdataPtr).ScanPtr=''
                    if ptr_valid((*(*pState).readdataPtr).BasePtr) then *(*(*pState).readdataPtr).BasePtr=''
                    if ptr_valid((*(*pState).readdataPtr).TypePtr) then *(*(*pState).readdataPtr).TypePtr=''
                    if ptr_valid((*(*pState).readdataPtr).Collimationptr) then *(*(*pState).readdataPtr).Collimationptr=''
                    if ptr_valid((*(*pState).readdataPtr).MosaicPtr) then *(*(*pState).readdataPtr).MosaicPtr=''
                    if ptr_valid((*(*pState).readdataPtr).XtalOrientPtr) then *(*(*pState).readdataPtr).XtalOrientPtr=''
                    if ptr_valid((*(*pState).readdataPtr).LatticeAPtr) then *(*(*pState).readdataPtr).LatticeAPtr=''
                    if ptr_valid((*(*pState).readdataPtr).LatticeBPtr) then *(*(*pState).readdataPtr).LatticeBPtr=''
                    if ptr_valid((*(*pState).readdataPtr).LatticeCPtr) then *(*(*pState).readdataPtr).LatticeCPtr=''
                    if ptr_valid((*(*pState).readdataPtr).LatticeAlPtr) then *(*(*pState).readdataPtr).LatticeAlPtr=''
                   if ptr_valid( (*(*pState).readdataPtr).LatticeBePtr) then *(*(*pState).readdataPtr).LatticeBePtr=''
                    if ptr_valid((*(*pState).readdataPtr).LatticeGaPtr) then *(*(*pState).readdataPtr).LatticeGaPtr=''
                    if ptr_valid((*(*pState).readdataPtr).E_ZEROPtr) then *(*(*pState).readdataPtr).E_ZEROPtr=''
                   if ptr_valid( (*(*pState).readdataPtr).DELTA_EPtr) then *(*(*pState).readdataPtr).DELTA_EPtr=''
                    if ptr_valid((*(*pState).readdataPtr).E_FINALPtr) then *(*(*pState).readdataPtr).E_FINALPtr=''
                    if ptr_valid((*(*pState).readdataPtr).MON_D_SPACEPtr) then *(*(*pState).readdataPtr).MON_D_SPACEPtr=''
                    if ptr_valid((*(*pState).readdataPtr).AN_D_SPACEPtr) then *(*(*pState).readdataPtr).AN_D_SPACEPtr=''
                    if ptr_valid((*(*pState).readdataPtr).T_startPtr) then *(*(*pState).readdataPtr).T_startPtr=''
                    if ptr_valid((*(*pState).readdataPtr).T_INCPtr) then *(*(*pState).readdataPtr).T_INCPtr=''
                    if ptr_valid((*(*pState).readdataPtr).HH_Ptr) then *(*(*pState).readdataPtr).HH_Ptr=''
                    if ptr_valid((*(*pState).readdataPtr).KK_Ptr) then *(*(*pState).readdataPtr).KK_Ptr=''
                    if ptr_valid((*(*pState).readdataPtr).LL_Ptr) then *(*(*pState).readdataPtr).LL_Ptr=''
                    if ptr_valid((*(*pState).readdataPtr).dH_Ptr) then *(*(*pState).readdataPtr).dH_Ptr=''
                    if ptr_valid((*(*pState).readdataPtr).dK_Ptr) then *(*(*pState).readdataPtr).dK_Ptr=''
                    if ptr_valid((*(*pState).readdataPtr).dL_Ptr) then *(*(*pState).readdataPtr).dL_Ptr=''
                    if ptr_valid((*(*pState).readdataPtr).HFIELD_Ptr) then *(*(*pState).readdataPtr).HFIELD_Ptr=''
                    if ptr_valid((*(*pState).readdataPtr).MonitorPtr) then *(*(*pState).readdataPtr).MonitorPtr=''
                    if ptr_valid((*(*pState).readdataPtr).PrefactorPtr) then *(*(*pState).readdataPtr).PrefactorPtr=''
                    if ptr_valid((*(*pState).readdataPtr).QxPtr) then *(*(*pState).readdataPtr).QxPtr=''
                    if ptr_valid((*(*pState).readdataPtr).QyPtr) then *(*(*pState).readdataPtr).QyPtr=''
                    if ptr_valid((*(*pState).readdataPtr).QzPtr) then *(*(*pState).readdataPtr).QzPtr=''
                    if ptr_valid((*(*pState).readdataPtr).EnergyPtr) then *(*(*pState).readdataPtr).EnergyPtr=''
                    if ptr_valid((*(*pState).readdataPtr).TempPtr) then *(*(*pState).readdataPtr).TempPtr=''
                    if ptr_valid((*(*pState).readdataPtr).TimePtr) then *(*(*pState).readdataPtr).TimePtr=''
                    if ptr_valid((*(*pState).readdataPtr).TotalIntPtr) then *(*(*pState).readdataPtr).TotalIntPtr=''

 widget_control, event.top, set_uvalue = pState

file = (*pstate).dfile

 f2=StrSplit(file,'/',/Extract)
 file2=f2[N_elements(f2)-1]
 wdave2=where(strmatch(file2,'*.dave',/fold_case) eq 1,wheredave)
 wbt4=where(strmatch(file2,'*.bt4',/fold_case) eq 1,wherebt4)
 wbta=where(strmatch(file2,'*.bta',/fold_case) eq 1,wherebta)

if not wheredave and not wherebt4 and not wherebta then return

if wherebta then gaussReadbta,event
if wherebt4 then gaussReadbt4,event
if wheredave then gaussReaddave,event

endif
end
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gaussplotcalc,event
    widget_control,event.top,get_uvalue = pState
    if (ptr_valid(pstate) eq 0) then return
wset,(*pState).winpix

;if (can show widget alive and get values) then begin
;
;endif


pfund=((*pState).plotfund)
pover=((*pState).plotover)
pdata=((*pState).plotdata)

treatment=['Filename: ' +(*pState).gfile]
treatment=[[treatment],'Gaussian job run: '+strjoin((*pState).job)]
treatment=[[treatment],'max. force= '+string(((*pState).maxforce)[0])]
treatment=[[treatment],'rms. force= '+string(((*pState).rmsforce)[0])]
treatment=[[treatment],'max. disp= '+string(((*pState).maxdisp)[0])]
treatment=[[treatment],'rms. disp= '+string(((*pState).rmsdisp)[0])]
treatment=[[treatment],'numberatoms= '+string((*pState).numberatoms)]

; define a circular symbol for plotting!!!! puts in psym=8
th=hfbs_makepoints(xlo=0.0,xhi=2.0*!pi,npts=20)
xsym=cos(th)& ysym=sin(th)
usersym,xsym,ysym


;check the delta function here:
mono=widget_info((*pState).monof, /DROPLIST_SELECT)

if mono eq 2 then begin
    x=(*(*pState).origfreqptr)
    widget_control,(*pstate).mevbutton,get_value=mevflag
    if (mevflag ne 'Freq. in wavenumber') then begin
         x=x/8.0655
    endif
    y=*(*pstate).eigenyptr
    xmin=(*pstate).xminrange
    xmax=(*pstate).xmaxrange
    ymin=(*pstate).yminrange
    ymax=(*pstate).ymaxrange
    freqstring='Energy Transfer (meV)'
    xunits='meV'
    treatment=[[treatment],'Spectra plotted as delta functions']

    if (mevflag eq 'Freq. in wavenumber') then begin
       freqstring='Energy Transfer (cm-1)'
       xunits='cm-1'
    endif
    plot, transpose([[x],[x],[x]]), transpose([[y*0],[y],[y*0]]),$
       xrange=[xmin,xmax],$;yrange=[ymin,ymax],
            ytitle='Intensity',xtitle=freqstring
widget_control,event.top,set_uvalue=pState,/no_copy
widget_control,event.top,get_uvalue = pState,/no_copy
wset,(*pState).winVis
device,copy = [0,0,!d.x_size,!d.y_size,0,0,(*pState).winPix]

if ((*pstate).dump eq 1) then begin
      file = DIALOG_PICKFILE(/WRITE, FILTER = '*.out')
   if n_elements(file) ne 0 then begin
        OPENW, datout, file, /GET_LUN
    FOR i=0,n_elements(x)-1 DO BEGIN
        printf,datout,x[i],y[i],format='(f10.4,x,f10.4)'
    ENDFOR
    FREE_LUN, datout

    OPENW, irout, file+'.ir', /GET_LUN
    FOR i=0,n_elements((*(*pState).origfreqptr))-1 DO BEGIN
        printf,irout,(*(*pState).origfreqptr)[i],(*(*pState).iriptr)[i],format='(f10.4,x,f10.4)'
    ENDFOR
    FREE_LUN, irout
   endif
(*pstate).dump=0
endif

if ((*pstate).dumpdave eq 1) then begin
      fileout = DIALOG_PICKFILE(/WRITE, FILTER = '*.dave')
   if (fileout ne '') then begin

f2=StrSplit(fileout,'/',/Extract)
file2=f2[N_elements(f2)-1]
wheredave2=where(strmatch(file2,'*.dave',/fold_case) eq 1,dave2)
if dave2 ge 1 then begin
   file=fileout
endif else begin
  file=fileout+'_g98.dave'
endelse

    ret=create_dave_pointer(daveptr,$
          INSTRUMENT = 'GAUSSIAN',$
          QTY = y,$
          QTUNITS = 'Calculated Intensity:arb_units',$
          QTLABEL = 'I(Calc) (arb units)',$
          XVALS = x,$
          XTYPE = 'POINTS',$
          XUNITS = xunits,$
          XLABEL = 'Energy Transfer',$
          ;SPECIFICSTR = specificstr,$
          TREATMENT = treatment,$
          ERMSG = errmsg)
    if ret eq 0 then begin
       msg = 'Problem creating the DAVE pointer'
       if ptr_valid(daveptr) then heap_free,daveptr
       (*pstate).dumpdave=0
       widget_control,event.top,set_uvalue = pState,/no_copy
       return
    endif
    if ptr_valid(daveptr) then begin
       save,daveptr,filename = file
            heap_free,daveptr
    endif
     endif
(*pstate).dumpdave=0
endif

widget_control,event.top,set_uvalue = pState,/no_copy
return
end


if (pfund eq 1) then begin
    y=((*pState).fundamental)
    treatment=[[treatment],'Calculated spectra from Fundamentals']
endif
if (pover eq 1) then begin
    y=((*pState).overtone)
    treatment=[[treatment],'Calculated spectrum has Overtones included']
endif
if (pdata eq 1) then begin
   if (n_elements(*(*(*pState).readdataPtr).Energyptr) ne 0) then xdata=*(*(*pState).readdataPtr).EnergyPtr
   if (n_elements(*(*(*pState).readdataPtr).TotalIntPtr) ne 0) then ydata=*(*(*pState).readdataPtr).TotalIntPtr
endif else y=(*pState).fundamental

if ((pover eq 1) and (pfund eq 1)) then begin
    y=((*pState).overtone)+ (*pState).fundamental

endif

if ((pover eq 1) OR (pfund eq 1)) then begin
    x=((*pState).freqspec)

    xmin=(*pstate).xminrange
    xmax=(*pstate).xmaxrange
    ymin=(*pstate).yminrange
    ymax=(*pstate).ymaxrange

    freqstring='Energy Transfer (meV)'
    widget_control,(*pstate).mevbutton,get_value=mevflag
    if (mevflag eq 'Freq. in wavenumber') then begin
       freqstring='Energy Transfer (cm-1)'
       xunits='cm-1'
    endif
    plot,x,y,xrange=[xmin,xmax],yrange=[ymin,ymax],xstyle=1,ystyle=1,$
       ytitle='Intensity',xtitle=freqstring;,psym=calcsymval

    if (mevflag eq 'Freq. in wavenumber') then begin
        if (pdata eq 1) and (n_elements(xdata) ge 2) then begin
            xxdata=xdata*8.06544
           oplot,xxdata,ydata,psym=8,symsize=.7
        endif
   endif  else begin
        if  (pdata eq 1) and (n_elements(xdata) ge 1) then begin
            oplot,(xdata),ydata,psym=8,symsize=.7


         endif
   endelse
endif else begin
        if  (pdata eq 1) and (n_elements(xdata) ge 1) then begin
            plot,(xdata),ydata,psym=datasymval,symsize=.7

         endif
endelse


widget_control,event.top,set_uvalue=pState,/no_copy
widget_control,event.top,get_uvalue = pState,/no_copy
wset,(*pState).winVis
device,copy = [0,0,!d.x_size,!d.y_size,0,0,(*pState).winPix]
if ((*pstate).dump eq 1) then begin
      file = DIALOG_PICKFILE(/WRITE, FILTER = '*.out')
   if (file ne '') then begin
        OPENW, datout, file, /GET_LUN
    FOR i=0,n_elements(x)-1 DO BEGIN
        printf,datout,x[i],y[i],format='(f10.4,x,f10.4)'
    ENDFOR
    FREE_LUN, datout

    OPENW, irout, file+'.ir', /GET_LUN
    FOR i=0,n_elements((*(*pState).origfreqptr))-1 DO BEGIN
        printf,irout,(*(*pState).origfreqptr)[i],(*(*pState).iriptr)[i],format='(f10.4,x,f10.4)'
    ENDFOR
    FREE_LUN, irout
   endif
(*pstate).dump=0
endif
;;;;;;;;;; add a dave file dump option.........

if ((*pstate).dumpdave eq 1) then begin
      fileout = DIALOG_PICKFILE(/WRITE, FILTER = '*.dave')
   if (fileout ne '') then begin

f2=StrSplit(fileout,'/',/Extract)
file2=f2[N_elements(f2)-1]
wheredave2=where(strmatch(file2,'*.dave',/fold_case) eq 1,dave2)
if dave2 ge 1 then begin
   file=fileout
endif else begin
  file=fileout+'_g98.dave'
endelse


    ret=CREATE_DAVE_POINTER(daveptr,$
          INSTRUMENT = 'GAUSSIAN',          $
          QTY = y,       $
          QTUNITS = 'Calculated Intensity:arb_units', $
          QTLABEL = 'I(Calc) (arb units)',         $
          XVALS = x,             $
          XTYPE = 'POINTS',             $
          XUNITS = xunits,              $
          XLABEL = 'Energy Transfer',           $
          ;SPECIFICSTR = specificstr,        $
          TREATMENT = treatment,   $
          ERMSG = errmsg             )
    if ret eq 0 then begin
       msg = 'Problem creating the DAVE pointer'
       if ptr_valid(daveptr) then heap_free,daveptr
       (*pstate).dumpdave=0
       widget_control,event.top,set_uvalue = pState,/no_copy
       return
    endif
    if ptr_valid(daveptr) then begin
       save,daveptr,filename = file
            heap_free,daveptr
    endif
     endif
(*pstate).dumpdave=0
endif

widget_control,event.top,set_uvalue = pState,/no_copy
return
end
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gaussReadGaussian,event

      widget_control,event.top,get_uvalue = pState

rotconstants=(*pState).rotconstants
lowfrequencies=(*pState).lowfrequencies
job=((*pState).job)
dof=(*pState).dof
file=(*pState).gfile


maxforce=(*pState).maxforce
maxforceflag=(*pState).maxforceflag
rmsforce=(*pState).rmsforce
rmsforceflag=(*pState).rmsforceflag
maxdisp=(*pState).maxdisp
maxdispflag=(*pState).maxdispflag
rmsdisp=(*pState).rmsdisp
rmsdispflag=(*pState).rmsdispflag
numberatoms=(*pState).numberatoms

    nbig=50000l        ;max size file and text array
      dumstr=''          ;dummy read in variable
      Numlines=0         ;initialize the total number of lines
      startline=0      ;initialize the line we want to read from
      jobline=0          ;this is to extract the G98 job types
      dofline=0          ;this is for the VIB. DEGREES OF FREEDOM
      rottline=0         ;this is for the ROTATIONAL TEMP(K)
; NB 11.6045 K =1meV and  1meV = 8.0655 cm-1, and 1 meV=0.2417 THz
      lowfreqline=0        ;this is for the low frequencies..
      convergeline=0         ;this is for the RMS factors..
      atomline=0         ;this is for the atom counting.
      ORIENTLINE=0       ;Input orientation counting

      textarr=strarr(nbig)   ;large text array

        openr,lun,file,/get_lun, error=err
        if (err ne 0) then print,'file isnt there'
        if (err ne 0) then return

        openr,lun,file,/get_lun, error=err
        tmp=0
        WHILE (NOT EOF(lun)) DO BEGIN
            readf,lun,dumstr
       textarr(tmp)=dumstr
            if (strcmp(strmid(dumstr,1,2),'# ')) then jobline=tmp
            if (strcmp(strmid(dumstr,1,15),'Deg. of freedom')) then dofline=tmp
            if (strcmp(strmid(dumstr,1,23),'ROTATIONAL TEMPERATURES'))$
            or (strcmp(strmid(dumstr,1,23),'Rotational temperatures'))$
           then rottline=tmp
            if (strcmp(strmid(dumstr,1,15),'Low frequencies')) then lowfreqline=tmp
       if (strcmp(strmid(dumstr,1,13),'Maximum Force')) $
         then convergeline=tmp
       if (strcmp(strmid(dumstr,1,20),'Total atomic charges')) $
          or (strcmp(strmid(dumstr,1,23),'Mulliken atomic charges')) $
         then atomline=tmp
        if (stregex(dumstr,$
                     'and normal coordinates',$
                      /boolean) eq 1) then startline=tmp
         if (stregex(dumstr,$
                     'Input orientation',$
                      /boolean) eq 1) then ORIENTLINE=tmp

       tmp=tmp+1
        ENDWHILE
        numlines=tmp
        free_lun,lun

;extract those things that may prove useful:
if (jobline ne 0) then begin
    if (n_elements(textarr(jobline)) ge 1) then $
       job=(strsplit(textarr(jobline),' ',/extract))
endif

if (dofline ne 0) then begin
    if (n_elements(textarr(dofline)) ge 1) then $
       dofa=(strsplit(textarr(dofline),' ',/extract))
       dof=fix(dofa[3])          ;convert deg. freedom to int.
endif

;dof=9    ; XXXXXXXXXXXXXXXXXXXXXXX


if (dofline ne 0) then begin
    if (n_elements(textarr(rottline)) ge 1) then $
       rota=(strsplit(textarr(rottline),' ',/extract))
    if (n_elements((rota)) ge 4) then $
       rotconstants[0]=rota(3)
    if (n_elements((rota)) ge 5) then $
       rotconstants[1]=rota(4)
    if (n_elements((rota)) ge 6) then $
       rotconstants[2]=rota(5)
endif

if (lowfreqline ne 0) then begin
    if (n_elements(textarr(lowfreqline-1)) ge 1) then $
       lfl=(strsplit(textarr(lowfreqline-1),' ',/extract))
    if (n_elements(lfl) ge 4) then $
       lowfrequencies[0]=lfl(3)
    if (n_elements(lfl) ge 5) then $
       lowfrequencies[1]=lfl(4)
    if (n_elements(lfl) ge 6) then $
       lowfrequencies[2]=lfl(5)
    if (n_elements(lfl) ge 7) then $
       lowfrequencies[3]=lfl(6)
    if (n_elements(lfl) ge 8) then $
       lowfrequencies[4]=lfl(7)
    if (n_elements(lfl) ge 9) then $
       lowfrequencies[5]=lfl(8)

    if (n_elements(textarr(lowfreqline)) ge 1) then $
       lfl=(strsplit(textarr(lowfreqline),' ',/extract))
    if (n_elements(lfl) ge 4) then $
       lowfrequencies[6]=lfl(3)
    if (n_elements(lfl) ge 5) then $
       lowfrequencies[7]=lfl(4)
    if (n_elements(lfl) ge 6) then $
       lowfrequencies[8]=lfl(5)
endif

if (convergeline ne 0) then begin
    maxf=(strsplit(textarr(convergeline),' ',/extract))
    if (n_elements(maxf) ge 2) then $
       maxforce[0]=maxf(2)
    if (n_elements(maxf) ge 3) then $
       maxforce[1]=maxf(3)
    if (n_elements(maxf) ge 4) then $
       maxforceflag=maxf(4)

    rmsf=(strsplit(textarr(convergeline+1),' ',/extract))
    if (n_elements(rmsf) ge 2) then $
       rmsforce[0]=rmsf(2)
    if (n_elements(rmsf) ge 3) then $
       rmsforce[1]=rmsf(3)
    if (n_elements(rmsf) ge 4) then $
       rmsforceflag=rmsf(4)

    maxd=(strsplit(textarr(convergeline+2),' ',/extract))
    if (n_elements(maxd) ge 2) then $
       maxdisp[0]=maxd(2)
    if (n_elements(maxd) ge 3) then $
       maxdisp[1]=maxd(3)
    if (n_elements(maxd) ge 4) then $
       maxdispflag=maxd(4)

    rmsd=(strsplit(textarr(convergeline+3),' ',/extract))
    if (n_elements(rmsd) ge 2) then $
       rmsdisp[0]=rmsd(2)
    if (n_elements(rmsd) ge 3) then $
       rmsdisp[1]=rmsd(3)
    if (n_elements(rmsd) ge 4) then $
       rmsdispflag=rmsd(4)


       if (STRMATCH(rmsdispflag , 'NO' )) then $
          print, 'OOPS!!! best check the RMS DISP. convergence criteria.'
       if (STRMATCH(maxdispflag , 'NO' )) then $
          print, 'OOPS!!! best check the MAX DISP. convergence criteria.'
       if (STRMATCH(rmsforceflag , 'NO' )) then $
          print, 'OOPS!!! best check the RMS FORCE convergence criteria.'
       if (STRMATCH(maxforceflag , 'NO' )) then $
          print, 'OOPS!!! best check the MAX FORCE convergence criteria.'

endif

if (n_elements(textarr(atomline)) ge 1) AND (atomline ne 0) then begin
       i=0
       while $
       ((strcmp(strmid(textarr(atomline+i+2),1,23),'Sum of Mulliken charges'))$
          ne 1) do begin
       tmp=(strsplit(textarr(atomline+i+2),' ',/extract))
       i=i+1
       endwhile
       numberatoms=i
       if (numberatoms eq 0) then print,'You have no atoms? Really?'
endif else begin
       print, 'This program seems to be confused. Where are the header lines'
       return
endelse


atomlist=strarr(numberatoms)
newatomlist=strarr(numberatoms)
atomicnumber=intarr(numberatoms)
newatomicnumber=intarr(numberatoms)
ixsection=fltarr(numberatoms)
newixsection=fltarr(numberatoms)
ixsection=fltarr(numberatoms)
mass=fltarr(numberatoms)
initialorientation=fltarr(numberatoms,3)

if (orientline ne 0) then begin
         i=0
         while $
         ((strcmp(strmid(textarr(orientline+i+5),1,23),'-----------------------'))$
          ne 1) do begin
         oline=(strsplit(textarr(orientline+i+5),' ',/extract))
         initialorientation[i,0]=oline[3]
         initialorientation[i,1]=oline[4]
         initialorientation[i,2]=oline[5]
         i=i+1
         endwhile

endif

if (n_elements(textarr(atomline)) ge 1) then begin
       i=0
       while $
       ((strcmp(strmid(textarr(atomline+i+2),1,23),'Sum of Mulliken charges'))$
          ne 1) do begin
       tmp=(strsplit(textarr(atomline+i+2),' ',/extract))
       atomlist[i]=tmp[1]
       i=i+1
       endwhile
       if (n_elements(atomlist) eq 0) then print,'You have no atoms? Really?'
endif


startline=startline+2
i=0
dof=1000
;make an arb number of jumps....

symmode=strarr(dof)
freq=fltarr(dof)
redmass=fltarr(dof)
frcc=fltarr(dof)
iri=fltarr(dof)
raman=fltarr(dof)
depolar=fltarr(dof)
vector=fltarr(dof,numberatoms,3)       ;put the displacements in x,y,z form

while (stregex(textarr(startline),'------',/boolean) eq 0) do begin

    symline=(strsplit(textarr(startline),' ',/extract))
    if (n_elements(symline) eq 3) then begin
       symmode[i]=symline[0]
           symmode[i+1]=symline[1]
       symmode[i+2]=symline[2]
       startline=startline+1
      if (strcmp(strmid(textarr(startline),1,14),'Frequencies --')) then begin
          freqline=(strsplit(textarr(startline),' ',/extract))
          freq[i]=freqline[2]
          freq[i+1]=freqline[3]
          freq[i+2]=freqline[4]
          startline=startline+1
       endif
       if (strcmp(strmid(textarr(startline),1,14),'Red. masses --')) then begin
          redline=(strsplit(textarr(startline),' ',/extract))
          redmass[i]=redline[3]
          redmass[i+1]=redline[4]
          redmass[i+2]=redline[5]
          startline=startline+1
       endif
       if (strcmp(strmid(textarr(startline),1,14),'Frc consts  --')) then begin
          frcline=(strsplit(textarr(startline),' ',/extract))
          frcc[i]=frcline[3]
          frcc[i+1]=frcline[4]
          frcc[i+2]=frcline[5]
          startline=startline+1
       endif
       if (strcmp(strmid(textarr(startline),1,14),'IR Inten    --')) then begin
          irline=(strsplit(textarr(startline),' ',/extract))
          irI[i]=irline[3]
          iri[i+1]=irline[4]
           iri[i+2]=irline[5]
          startline=startline+1
       endif
       if (strcmp(strmid(textarr(startline),1,14),'Raman Activ --')) then begin
          ramanline=(strsplit(textarr(startline),' ',/extract))
          raman[i]=ramanline[3]
          raman[i+1]=ramanline[4]
          raman[i+2]=ramanline[5]
          startline=startline+1
       endif
       if (strcmp(strmid(textarr(startline),1,14),'Depolar     --')) then begin
          depolarline=(strsplit(textarr(startline),' ',/extract))
          depolar[i]=depolarline[2]
          depolar[i+1]=depolarline[3]
          depolar[i+2]=depolarline[4]
          startline=startline+1
        endif
          startline=startline+1

       for j=0,numberatoms-1 do begin
         atdisp=(strsplit(textarr(startline),' ',/extract))
         atomicnumber[j]=atdisp[1]
         vector[i,j,0]=atdisp[2]
         vector[i,j,1]=atdisp[3]
         vector[i,j,2]=atdisp[4]
         vector[i+1,j,0]=atdisp[5]
         vector[i+1,j,1]=atdisp[6]
         vector[i+1,j,2]=atdisp[7]
         vector[i+2,j,0]=atdisp[8]
         vector[i+2,j,1]=atdisp[9]
         vector[i+2,j,2]=atdisp[10]
         startline=startline+1
       endfor

       startline=startline+1
       i=i+3
    endif else begin
        symmode[i]=symline[0]
           symmode[i+1]=symline[1]
       startline=startline+1
                freqline=(strsplit(textarr(startline),' ',/extract))
       freq[i]=freqline[2]
       freq[i+1]=freqline[3]
       startline=startline+1
                redline=(strsplit(textarr(startline),' ',/extract))
       redmass[i]=redline[3]
       redmass[i+1]=redline[4]
       startline=startline+1
                frcline=(strsplit(textarr(startline),' ',/extract))
       frcc[i]=frcline[3]
       frcc[i+1]=frcline[4]
       startline=startline+1
                irline=(strsplit(textarr(startline),' ',/extract))
       irI[i]=irline[3]
       iri[i+1]=irline[4]
       startline=startline+1
                ramanline=(strsplit(textarr(startline),' ',/extract))
       raman[i]=ramanline[3]
       raman[i+1]=ramanline[4]
       startline=startline+1
                depolarline=(strsplit(textarr(startline),' ',/extract))
       depolar[i]=depolarline[2]
       depolar[i+1]=depolarline[3]
       depolar[i+2]=depolarline[4]
       startline=startline+2
       for j=0,numberatoms-1 do begin
         atdisp=(strsplit(textarr(startline),' ',/extract))
         atomicnumber[j]=atdisp[1]
         vector[i,j,0]=atdisp[2]
         vector[i,j,1]=atdisp[3]
         vector[i+1,j,0]=atdisp[5]
         vector[i+1,j,1]=atdisp[6]
         vector[i+2,j,0]=atdisp[8]
         vector[i+2,j,1]=atdisp[9]
         startline=startline+1
       endfor
    ENDELSE

if (n_elements(freq) eq 0) then begin
    print,'Either this program is having problems'
    print,'.....or you dont have any freq info. Try again'
    return
endif
endwhile




isotopeinfo=(*(*pState).isotopeinfoptr)

for i=0,numberatoms-1 do begin
    j=0
    repeat begin
       test=strsplit((isotopeinfo[j]),',',/extract)
       j=j+1
    endrep until atomlist[i] eq test[0]

    mass[i]=test[1]
    ixsection[i]=test[4]
endfor

;how many elements?


max=n_elements(where(freq))-1
*(*pState).initialorientationptr=initialorientation
(*pState).rotconstants=rotconstants
(*pState).lowfrequencies=lowfrequencies
(*pState).job=job
(*pState).dof=max+1
(*pState).maxforce=maxforce
(*pState).maxforceflag=maxforceflag
(*pState).rmsforce=rmsforce
(*pState).rmsforceflag=rmsforceflag
(*pState).maxdisp=maxdisp
(*pState).maxdispflag=maxdispflag
(*pState).rmsdisp=rmsdisp
(*pState).rmsdispflag=rmsdispflag
(*pState).numberatoms=numberatoms
(*((*pState).atomorig).listptr)=atomlist
(*((*pState).atomnew).newlistptr)=atomlist
(*((*pState).atomorig).numberptr)=atomicnumber
(*((*pState).atomnew).newnumberptr)=atomicnumber
(*((*pState).atomorig).ixsectionptr)=ixsection
(*((*pState).atomnew).newixsectionptr)=ixsection
(*((*pState).atomorig).massptr)=mass
(*((*pState).atomnew).newmassptr)=mass
(*(*pState).symmodeptr)=symmode[0:max]
(*(*pState).origfreqptr)=freq[0:max]
(*(*pState).freqptr)=freq[0:max]
(*(*pState).redmassptr)=redmass[0:max]
(*(*pState).frccptr)=frcc[0:max]
(*(*pState).iriptr)=iri[0:max]
(*(*pState).ramanptr)=raman[0:max]
(*(*pState).depolarptr)=depolar[0:max]
(*(*pState).vectorptr)=vector[0:max,*,*]
      widget_control,event.top,set_uvalue=pState,/no_copy

gausseigenvector,event
end
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gaussReadbta,event
widget_control, event.top, get_uvalue = pState

  catch,errno
  if (errno ne 0) then begin
    catch,/cancel
    void = dialog_message(/error,["Problem with your BTA file?",!err_string])
    return
  endif

dfile=(*pState).dfile
                 linenum=0
                  ndets = 128             ;number detectors
                  str = strarr(12)        ;keep header info here
                  dumstr=''               ;dummy
                  nbig=5000               ;More lines to read than neccessary
                  DATE=''                 ;The file Date
                  SCAN=''                 ;The scan type.
                  monitor=0.0             ;monitor counts
                  prefactor=0.0           ;and how many prefactors look at a BT4 file for the rest...
                  base=''
                  type=''
                  Collimation=''
                  Mosaic=''
                  XtalOrient=''
                  LatticeA=''
                  LatticeB=''
                  LatticeC=''
                  LatticeAl=''
                  LatticeBe=''
                  LatticeGa=''
                  E_ZERO=0.0
                  DELTA_E=0.0
                  E_FINAL=0.0
                  MON_D_SPACE=0.0
         AN_D_SPACE=0.0
                  T_Start=0.0
                  T_Inc=0.0
                  H=0.0
                  K=0.0
                  L=0.0
                  dH=0.0
                  dK=0.0
                  dL=0.0
                  HFIELD=0.0
                  T_flag = 0              ;is there Temperature recorded?
                  Qx=fltarr(nbig)
                  Qy=fltarr(nbig)
                  Qz=fltarr(nbig)
                  E=fltarr(nbig)
                  TEMP=fltarr(nbig)
                  TIME=fltarr(nbig)
                  INTENSITY=fltarr(nbig)
                  ;##########################

openr,lun,dfile,/get_lun
         ; Read the header
                  for j = 0,11 do begin
                  readf,lun,dumstr
                  str[linenum]=dumstr
                  linenum=linenum+1
                  endfor
                  ;extract everything incase I want it later.
                  Date=(STRMID(str[0],16,17))
                  scan=(STRMID(str[0],36,1))
                  monitor=float(STRMID(str[0],39,12))
                  prefactor=float(STRMID(str[0],52,4))
                  base=(STRMID(str[0],58,4))
                  type=(STRMID(str[0],71,3))
                    Collimation=(STRMID(str[3],1,11))
                    Mosaic=(STRMID(str[3],18,8))
                    XtalOrient=(STRMID(str[3],36,23))
                    LatticeA=float(strtrim((STRMID(str[5],1,9)),2))
                    LatticeB=float(strtrim((STRMID(str[5],11,9)),2))
                    LatticeC=float(strtrim((STRMID(str[5],21,9)),2))
                    LatticeAl=float(strtrim((STRMID(str[5],33,7)),2))
                    LatticeBe=float(strtrim((STRMID(str[5],43,7)),2))
                    LatticeGa=float(strtrim((STRMID(str[5],53,7)),2))
                    E_ZERO=float(strtrim((STRMID(str[7],3,7)),2))
                    DELTA_E=float(strtrim((STRMID(str[7],13,7)),2))
                    E_FINAL=float(strtrim((STRMID(str[7],23,10)),2))
                    MON_D_SPACE=float(strtrim((STRMID(str[7],41,5)),2))
                    AN_D_SPACE=float(strtrim((STRMID(str[7],49,5)),2))
                    T_start=float(strtrim((STRMID(str[7],58,9)),2))
                    T_INC=float(strtrim((STRMID(str[7],68,9)),2))
                    HH=float(strtrim((STRMID(str[9],0,7)),2))
                    KK=float(strtrim((STRMID(str[9],8,7)),2))
                    LL=float(strtrim((STRMID(str[9],16,7)),2))
                    dH=float(strtrim((STRMID(str[9],26,7)),2))
                    dK=float(strtrim((STRMID(str[9],34,7)),2))
                    dL=float(strtrim((STRMID(str[9],42,7)),2))
                    HFIELD=float(strtrim((STRMID(str[9],51,7)),2))

                  ;####################################


                  ;####################################
                  ;check for the extra column
                  if STRMID(str[11],37,1) eq 'T' then T_flag = 1
                  if T_flag eq 1 then begin
           line=0
                          count=0
                          WHILE (NOT EOF(lun)) DO BEGIN
                            READF, lun,dumstr
;add a clause to remove any old file fake lines with just a single zero in.
                                if (dumstr eq '  0') or (dumstr eq ' 0') $
                 or (dumstr eq '0') then begin
             line=line+1
             continue
          endif
             Qx(count)=float(STRMID(dumstr,0,8))
                            Qy(count)=float(STRMID(dumstr,8,8))
                            Qz(count)=float(STRMID(dumstr,16,8))
                            E(count)=float(STRMID(dumstr,24,9))
                            TEMP(count)=float(STRMID(dumstr,33,9))
                            TIME(count)=float(STRMID(dumstr,42,8))
                            INTENSITY(count)=float(STRMID(dumstr,53,11))
                            count=count+1
             line=line+1
                          ENDWHILE
                  ;###########################################
                  endif else begin
                         count=0
          line=0
                          WHILE (NOT EOF(lun)) DO BEGIN
                            READF, lun,dumstr
;add a clause to remove any old file fake lines with just a single zero in.
          if (dumstr eq '  0') or (dumstr eq ' 0') $
                            or (dumstr eq '0') then begin
             line=line+1
             continue
          endif
                            Qx(count)=float(STRMID(dumstr,0,8))
                            Qy(count)=float(STRMID(dumstr,8,8))
                            Qz(count)=float(STRMID(dumstr,16,8))
                            E(count)=float(STRMID(dumstr,24,9))
                            TEMP(count)=0.0
                            TIME(count)=float(STRMID(dumstr,33,8))
                            INTENSITY(count)=float(STRMID(dumstr,42,12))
                            count=count+1
             line=line+1

                         ENDWHILE
                  endelse

free_lun,lun

                  *(*(*pState).readdataPtr).DatePtr=date
                  *(*(*pState).readdataPtr).ScanPtr=scan
                  *(*(*pState).readdataPtr).BasePtr=base
                  *(*(*pState).readdataPtr).TypePtr=type
                  *(*(*pState).readdataPtr).Collimationptr=collimation
                  *(*(*pState).readdataPtr).MosaicPtr=mosaic
                  *(*(*pState).readdataPtr).XtalOrientPtr=xtalorient
                  *(*(*pState).readdataPtr).LatticeAPtr=latticeA
                  *(*(*pState).readdataPtr).LatticeBPtr=latticeB
                  *(*(*pState).readdataPtr).LatticeCPtr=latticeC
                  *(*(*pState).readdataPtr).LatticeAlPtr=latticeAL
                  *(*(*pState).readdataPtr).LatticeBePtr=latticeBe
                  *(*(*pState).readdataPtr).LatticeGaPtr=latticeGa
                  *(*(*pState).readdataPtr).E_ZEROPtr=E_zero
                  *(*(*pState).readdataPtr).DELTA_EPtr=Delta_e
                  *(*(*pState).readdataPtr).E_FINALPtr=E_final
                  *(*(*pState).readdataPtr).MON_D_SPACEPtr=MON_D_SPACE
                  *(*(*pState).readdataPtr).AN_D_SPACEPtr=AN_D_SPACE
                  *(*(*pState).readdataPtr).T_startPtr=T_start
                  *(*(*pState).readdataPtr).T_INCPtr=T_INC
                  *(*(*pState).readdataPtr).HH_Ptr=HH
                  *(*(*pState).readdataPtr).KK_Ptr=KK
                  *(*(*pState).readdataPtr).LL_Ptr=LL
                  *(*(*pState).readdataPtr).dH_Ptr=dH
                  *(*(*pState).readdataPtr).dK_Ptr=dK
                  *(*(*pState).readdataPtr).dL_Ptr=dL
                  *(*(*pState).readdataPtr).HFIELD_Ptr=HFIELD
                  *(*(*pState).readdataPtr).MonitorPtr=Monitor
                  *(*(*pState).readdataPtr).PrefactorPtr=Prefactor
                  *(*(*pState).readdataPtr).QxPtr=qx[0:count-1]
                  *(*(*pState).readdataPtr).QyPtr=qy[0:count-1]
                  *(*(*pState).readdataPtr).QzPtr=qz[0:count-1]
                  *(*(*pState).readdataPtr).EnergyPtr=E[0:count-1]
                  *(*(*pState).readdataPtr).TempPtr=TEMP[0:count-1]
                  *(*(*pState).readdataPtr).TimePtr=TIME[0:count-1]
                  *(*(*pState).readdataPtr).TotalIntPtr=INTENSITY[0:count-1]


widget_control,event.top,set_uvalue=pState,/no_copy
end

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gaussReadbt4,event
widget_control, event.top, get_uvalue = pState
  catch,errno
  if (errno ne 0) then begin
    catch,/cancel
    void = dialog_message(/error,["Problem with your BT4 file?",!err_string])
    return
  endif
dfile=(*pState).dfile
                 linenum=0
                  ndets = 128             ;number detectors
                  str = strarr(12)        ;keep header info here
                  dumstr=''               ;dummy
                  nbig=5000               ;More lines to read than neccessary
                  DATE=''                 ;The file Date
                  SCAN=''                 ;The scan type.
                  monitor=0.0             ;monitor counts
                  prefactor=0.0           ;and how many prefactors look at a BT4 file for the rest...
                  base=''
                  type=''
                  Collimation=''
                  Mosaic=''
                  XtalOrient=''
                  LatticeA=''
                  LatticeB=''
                  LatticeC=''
                  LatticeAl=''
                  LatticeBe=''
                  LatticeGa=''
                  E_ZERO=0.0
                  DELTA_E=0.0
                  E_FINAL=0.0
                  MON_D_SPACE=0.0
         AN_D_SPACE=0.0
                  T_Start=0.0
                  T_Inc=0.0
                  H=0.0
                  K=0.0
                  L=0.0
                  dH=0.0
                  dK=0.0
                  dL=0.0
                  HFIELD=0.0
                  T_flag = 0              ;is there Temperature recorded?
                  Qx=fltarr(nbig)
                  Qy=fltarr(nbig)
                  Qz=fltarr(nbig)
                  E=fltarr(nbig)
                  TEMP=fltarr(nbig)
                  TIME=fltarr(nbig)
                  INTENSITY=fltarr(nbig)
                  ;##########################

openr,lun,dfile,/get_lun
        line = ''
    header = strarr(5)
    readf, lun, header
tflag=0
*(*(*pState).readdataPtr).DatePtr = strmid(header[0],16,17)
*(*(*pState).readdataPtr).ScanPtr = strmid(header[0],36,1)
*(*(*pState).readdataPtr).MonitorPtr = long(strmid(header[0],39,12))
*(*(*pState).readdataPtr).PrefactorPtr = long(strmid(header[0],52,4))
*(*(*pState).readdataPtr).BasePtr = strmid(header[0],58,4)
*(*(*pState).readdataPtr).typePtr= strmid(header[0],71,3)
*(*(*pState).readdataPtr).Collimationptr = strmid(header[3],1,11)
*(*(*pState).readdataPtr).MosaicPtr = strmid(header[3],18,8)
*(*(*pState).readdataPtr).XtalOrientPtr = strmid(header[3],36,23)
buf = fltarr(6) & readf,lun,buf
*(*(*pState).readdataPtr).LatticeAPtr = buf[0]
*(*(*pState).readdataPtr).LatticeBPtr = buf[1]
*(*(*pState).readdataPtr).LatticeCPtr = buf[2]
*(*(*pState).readdataPtr).LatticeAlPtr = buf[3]
*(*(*pState).readdataPtr).LatticeBePtr = buf[4]
*(*(*pState).readdataPtr).LatticeGaPtr = buf[5]
readf, lun,line & buf = fltarr(7) & readf, lun,buf
*(*(*pState).readdataPtr).E_ZEROPtr = buf[0]
*(*(*pState).readdataPtr).DELTA_EPtr = buf[1]
*(*(*pState).readdataPtr).E_FINALPtr = buf[2]
*(*(*pState).readdataPtr).MON_D_SPACEPtr = buf[3]
*(*(*pState).readdataPtr).AN_D_SPACEPtr = buf[4]
*(*(*pState).readdataPtr).T_startPtr = buf[5]
*(*(*pState).readdataPtr).T_INCPtr = buf[6]
readf, lun,line & readf, lun,buf
*(*(*pState).readdataPtr).HH_Ptr = buf[0]
*(*(*pState).readdataPtr).KK_Ptr = buf[1]
*(*(*pState).readdataPtr).LL_Ptr = buf[2]
*(*(*pState).readdataPtr).dH_Ptr = buf[3]
*(*(*pState).readdataPtr).dK_Ptr = buf[4]
*(*(*pState).readdataPtr).dl_Ptr = buf[5]
*(*(*pState).readdataPtr).HFIELD_Ptr = buf[6]
readf, lun,line,line
if (strmid(line,37,1) eq 'T') then tflag = 1


; now in the data section.
hbuf = (tflag)? fltarr(7) : fltarr(6)
dbuf = lonarr(ndets)
nspec = 2000 & i = 0
time = (temp = (totErrs = (totCnts =  (omega = (qx = (qy = (qz = fltarr(nspec))))))))
nosActDets = ndets
detCnts = fltarr(nosActDets,nspec)
detErrs = fltarr(nosActDets,nspec)

while (not eof(lun)) do begin
    readf, lun, hbuf,dbuf       ; read one data block (per energy transfer)
    qx[i] = hbuf[0]
    qy[i] = hbuf[1]
    qz[i] = hbuf[2]
    omega[i] = hbuf[3]
    if ( tflag) then begin
        temp[i] = hbuf[4]
        time[i] = hbuf[5]
        totCnts[i] = hbuf[6]
    endif else begin
        temp[i] = temp[i]
        time[i] = hbuf[4]
        totCnts[i] = hbuf[5]
    endelse
    totErrs[i] = (totCnts[i] le 0)? 1 : sqrt(totCnts[i]) ; assign 1 to erros if cnts <= 0
    detCnts[*,i] = dbuf[0:nosActDets-1] ; trim to active detectors only
    res = where(dbuf le 0)      ; locate indices where detcts <= 0
    if (res[0] ne -1) then begin ; and assign an error value of 1 to
        dbuf[res] = 1
    endif
    detErrs[*,i] = sqrt(dbuf[0:nosActDets-1])
    i = i + 1
endwhile
free_lun,lun

void=where(omega,count)

                  *(*(*pState).readdataPtr).QxPtr=qx[0:count-1]
                  *(*(*pState).readdataPtr).QyPtr=qy[0:count-1]
                  *(*(*pState).readdataPtr).QzPtr=qz[0:count-1]
                  *(*(*pState).readdataPtr).EnergyPtr=omega[0:count-1]
                  *(*(*pState).readdataPtr).TempPtr=TEMP[0:count-1]
                  *(*(*pState).readdataPtr).TimePtr=TIME[0:count-1]
                  *(*(*pState).readdataPtr).TotalIntPtr=totcnts[0:count-1]


widget_control,event.top,set_uvalue=pState,/no_copy
end

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gaussReaddave,event

  catch,errno
  if (errno ne 0) then begin
    catch,/cancel
    void = dialog_message(/error,["Problem with your DAVE file?",!err_string])
    return
  endif


widget_control, event.top, get_uvalue = pState
dfile=(*pState).dfile

    commonstr={instrument:'',$
                  histptr:ptr_new(/allocate_heap) , $
                  xlabel:'',$
                  ylabel:'',$
                  xtype:'',$
                  ytype:'',$
                  xunits:'',$
                  yunits:'',$
                  histlabel:'',$
                  histunits:'',$
                  treatmentptr:ptr_new(/allocate_heap) , $
                  commonstrptr:ptr_new(/allocate_heap)}



         datastr= {commonstr:commonstr, specificptr:ptr_new(specificstr,/no_copy)}
          datastrptr=ptr_new(datastr,/no_copy)

          descriptrstr={descriptrstr,$
                          name:'',$
                          units:'',$
                          legend:'',$
                          qty:0.0,$
                          err:0.0}

          davestr={daveptrstr,$
                  datastrptr:datastrptr, descriptr:ptr_new(descriptrstr)}
          davePtr=ptr_new(davestr,/no_copy)


restore,dfile

      y=(*(*(*davePtr).dataStrPtr).commonStr.histPtr).qty
      e=(*(*(*davePtr).dataStrPtr).commonStr.histPtr).err
      x=(*(*(*davePtr).dataStrPtr).commonStr.histPtr).x
      if (size(y))[0] eq 2 then begin
          y=y[*,0]
          e=e[*,0]
      endif

*(*(*pState).readdataPtr).EnergyPtr=x
*(*(*pState).readdataPtr).TotalIntPtr=y

    ptr_free,(*(*davePtr).datastrptr).commonstr.histptr
     ptr_free,(*(*davePtr).datastrptr).commonstr.commonstrptr
     ptr_free,(*(*davePtr).datastrptr).commonstr.treatmentptr
     ptr_free,(*davePtr).datastrptr
     ptr_free,(*davePtr).descriptr
     ptr_free,davePtr



widget_control,event.top,set_uvalue=pState,/no_copy
end

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gaussautoscaleplot,event
    widget_control,event.top,get_uvalue=pState
    if (n_elements(*(*pstate).freqptr) eq 0) then return
    xmin=(*pstate).xminrange
    xmax=(*pstate).xmaxrange
    ymin=(*pstate).yminrange
    ymax=(*pstate).ymaxrange
    pover=(*pstate).plotover
    pfund=(*pstate).plotfund
    pdata=(*pstate).plotdata

    freq=((*pState).freqspec)
    overtone=((*pState).overtone)
    fundamental=((*pState).fundamental)
    int=((*pState).plotspec)

xmin=min(freq)-10.0
xmax=max(freq)+10.0
if (pfund eq 1) then begin
    ymin=min(fundamental)-(2.0)
    ymax=max(fundamental)+(2.0)
endif
if (pover eq 1) then begin
    ymin=min(overtone)-(2.0)
    ymax=max(overtone)+(2.0)
endif

    if ((pfund eq 1) and (pover eq 1)) then begin
       if (ymin lt min(overtone+fundamental)) then begin
         ymin=ymin
       endif else begin
         ymin=min(overtone+fundamental)-2
       endelse
       if (ymax gt max(overtone+fundamental)+2) then begin
         ymax=ymax
       endif else begin
         ymax=max(overtone+fundamental)+2
       endelse
    endif

    (*pstate).xminrange=xmin
    (*pstate).xmaxrange=xmax
    (*pstate).yminrange=ymin
    (*pstate).ymaxrange=ymax
if (widget_info((*pState).tool2,/valid_id) ne 0) then begin
    if ((*pState).tool2 gt 0) then widget_control,(*pState).tool2,set_value=xmin
    if ((*pState).tool3 gt 0) then widget_control,(*pState).tool3,set_value=xmax
    if ((*pState).tool4 gt 0) then widget_control,(*pState).tool4,set_value=ymin
    if ((*pState).tool5 gt 0) then widget_control,(*pState).tool5,set_value=ymax
endif

    widget_control,event.top,set_uvalue=pState,/no_copy
    gaussplotcalc,event
end
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gaussmono_events,event
gaussresolution,event
end
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gausscalccolor_events,event

end
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gaussdatacolor_events,event

end

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

pro gausscoll_events,event
gaussresolution,event
end
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

pro gaussfinalspectrum,event
    widget_control,event.top,get_uvalue=pState
    spectrum=(*(*pState).spectrumptr)
    eigeny=(*(*pState).eigenyptr)
    freq=(*(*pState).freqptr)
    freqfactor=((*pState).freqfactor)
    freq=freq*freqfactor
    scalefactor=(*pstate).scalefactor
    B0=(*pstate).flatbg
    B1=    (*pstate).linearbg
    resfn=(*pstate).resfn
    dwf=(*pstate).dwf
    oversf=(*pstate).oversf

    resfn=resfn

    fundamental=(*pstate).fundamental
    overtone=(*pstate).overtone
    plotspec=(*pstate).plotspec
    freqspec=(*pstate).freqspec

    bg=freqspec
;need to get discrete resolutions at each freq
res=fltarr(n_elements(freq))

widget_control,(*pstate).mevbutton,get_value=mevflag

if (mevflag eq 'Freq. in wavenumber') then begin
  for i=0,n_elements(freq)-1 do begin
    k=0
    P = WHERE((resfn[*,0]) GE freq[i])
    if (mevflag eq 'Freq. in wavenumber') then begin
        if freq[i] le 2400 then begin
       res[i] =resfn[P[0],1]
    endif else begin
       res[i] =resfn[380,1]
    endelse
    endif else begin
        if freq[i] le 300 then begin
       res[i] =resfn[P[0],1]
    endif else begin
       res[i] =resfn[380,1]
    endelse
    endelse
  endfor
endif else begin
resfn=resfn/8.0655
 for i=0,n_elements(freq)-1 do begin
    k=0
    P = WHERE((resfn[*,0]) GE freq[i])
    if (mevflag ne 'Freq. in wavenumber') then begin
        if freq[i] le 297.56 then begin
       res[i] =resfn[P[0],1]
    endif else begin
       res[i] =resfn[380,1]
    endelse
    endif else begin
        if freq[i] le 37.19 then begin
       res[i] =resfn[P[0],1]
    endif else begin
       res[i] =resfn[380,1]
    endelse
    endelse
  endfor

endelse

;time0 = systime(1)
;for i=0,n1-1 do begin
;   tempI=0
;   for j=0,n2-1 do begin
;     TempI = TempI + (eigeny(j) / ((freqspec(i) - Freq(j))^2 + res(j)^2))
;   endfor
;
;   bg= 0.0 ; B0 + B1*i
;fundamental(i) = TempI*(exp(-1*dwf*dwf*freqspec(i)))*scalefactor+bg
;endfor
;print,'Loop way: ',(systime(1)-time0),' seconds'
;fundamental1=fltarr(n1)
;sum=fltarr(n1)
;tstart = systime(1)


n1 = n_elements(freqspec)
n2 = n_elements(eigeny)
un1 = 1+bytarr(n1)
un2 = 1+bytarr(n2)
bgarr=freqspec
bgarr=freqspec*B1+B0
sumf = total((un1#eigeny)*exp(-(freqspec#un1-un1#Freq)^2/(un1#res)^2),2)
fundamental = scalefactor*sumf*exp(-(dwf^2)*freqspec)+bgarr
;print,'Vectorized way: ',systime(1)-tstart, ' seconds'

    (*pstate).fundamental=fundamental
;   (*pstate).freqspec=freqspec

if ((*pstate).plotover eq 1) then begin

    ; calculate all possible combinations as being
    ; the  sum of frequency and mulitply intensities.
    CombNum = (n_elements(freq)^2 + n_elements(freq))/2
    CombFreq=fltarr(CombNum)
    CombInty=fltarr(CombNum)
       i = 0
        For j = 0,n_elements(freq)-1 do begin
            CombFreq(i) = Freq(j)*2
            CombInty(i) = eigeny(j)^2
            i = i + 1
            For k = (j + 1),n_elements(freq)-1 do begin
                CombFreq(i) = Freq(j) + Freq(k)
                CombInty(i) = eigeny(j) * eigeny(k)
                i = i + 1
            endfor
        endfor


    ;need to get discrete resolutions at each freq

    rescomb=fltarr(n_elements(combfreq))
if (mevflag eq 'Freq. in wavenumber') then begin
    for i=0,n_elements(combfreq)-1 do begin
    k=0
        P = WHERE((resfn[*,0]) GE combfreq[i])
       if combfreq[i] le 2400 then begin
         rescomb[i] =resfn[P[0],1]
       endif else begin
         rescomb[i] =resfn[380,1]
       endelse
    endfor
endif else begin
    for i=0,n_elements(combfreq)-1 do begin
    k=0
        P = WHERE((resfn[*,0]) GE combfreq[i])
       if combfreq[i] le 297.56 then begin
         rescomb[i] =resfn[P[0],1]
       endif else begin
         rescomb[i] =resfn[380,1]
       endelse
    endfor
endelse
n1 = n_elements(freqspec)
n2 = n_elements(CombInty)

;time0=systime(1)
;   for i=0,n1-1 do begin
;     tempI=0.0
;     for j=0,n2-1 do begin
;      TempI = TempI + (CombInty(j) / ((i - combFreq(j))^2 + rescomb(j)^2))
;     endfor
;     bg=0.0;   B0 + B1*i
;   overtone(i) = TempI*(exp(-1*dwf*dwf*i))*oversf;+bg
;
;   endfor
;print,'Loop way: ',(systime(1)-time0),' seconds'


overtone1=fltarr(n1)
sum=fltarr(n1)
tstart = systime(1)
un1 = 1+bytarr(n1)
un2 = 1+bytarr(n2)
sum = total((un1#CombInty)*exp(-(freqspec#un1-un1#combFreq)^2/(un1#rescomb)^2),2)
overtone1 = oversf*sum*exp(-(dwf^2)*freqspec)
;print,'Vectorized way: ',systime(1)-tstart, ' seconds'



(*pstate).overtone=overtone1
endif


;(*pstate).plotspec=fundamental
    widget_control,event.top,set_uvalue=pState
gaussautoscaleplot,event
end
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

pro gausseigenvector,event
;from the current mass, vector displacement and eigenvalue
;need to get rid of mass-weighted coordinates and
;calculate the real neutron intensity
widget_control,event.top,get_uvalue=pState
   mass=(*((*pState).atomnew).newmassptr)
    ixs=(*((*pstate).atomnew).newixsectionptr)
    fullvector=(*(*pstate).vectorptr)
    eigvec=fltarr((*pState).dof,(*pState).numberatoms) ; ARRAY SIZE FROM READ
    evecs=fltarr((*pState).dof,(*pState).numberatoms,3)
    totty=fltarr((*pState).dof)
    i=0
    for i=0,(*pState).numberatoms-1 do begin
       for j=0,(*pState).dof-1 do begin
          eigvec(j,i)=sqrt(fullvector(j,i,0)^2+fullvector(j,i,1)^2+fullvector(j,i,2)^2)
          evecs(j,i,0)=fullvector(j,i,0)*mass(i)
          evecs(j,i,1)=fullvector(j,i,1)*mass(i)
          evecs(j,i,2)=fullvector(j,i,2)*mass(i)

       endfor
    endfor


; Renormalize those eigenvectors:
for k=0,(*pState).dof-1 do begin
    tot=0
    esum=0.0
    for j=0,(*pState).numberatoms-1 do begin
       for i=0,2 do begin
         esum=esum+evecs(k,j,i)^2
       endfor

    endfor
    esum=sqrt(esum)
    for j=0,(*pState).numberatoms-1 do begin
       for i=0,2 do begin
         evecs(k,j,i)=evecs(k,j,i)/esum
       endfor
       eigvec(k,j)=sqrt(evecs(k,j,0)^2+evecs(k,j,1)^2+evecs(k,j,2)^2)
       fin=(eigvec(k,j)^2)*ixs(j)*exp(-1.0/mass(j))/mass(j)
       tot=tot+fin
    endfor
    totty[k]=tot
endfor

*(*pstate).eigenyptr=totty
*(*pstate).spectrumptr=totty*((*pstate).scalefactor)

widget_control,event.top,set_uvalue=pState,/no_copy
gaussresolution,event
end
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
function gausseditatoms,group_leader=base,in,out,changed
gaussvaredit,out,in,changed,group=base,y_scroll_size=(n_elements(in.list))
return,out
return,changed
end

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gaussresolution,event
    widget_control,event.top,get_uvalue=pState
if (n_elements(*(*pstate).iriptr) eq 0) then return
    d_cu=(*pstate).d_cu
    d_pg=(*pstate).d_pg
    dmono=d_cu
    coll1=20.0
    coll2=20.0
    ; get monochromator and collimation
    mono=widget_info((*pState).monof, /DROPLIST_SELECT)
    coll=widget_info((*pState).collimf, /DROPLIST_SELECT)
case mono of
    0: begin
           dmono=d_cu
       end
    1: begin
       dmono=d_pg
    end
    2: begin
       gaussplotcalc,event
       return
    end
    endcase

case coll of
    0: begin
           coll1=10.0
           coll2=10.0
       end
    1: begin
        coll1=10.0
        coll2=20.0
       end
    2: begin
        coll1=20.0
        coll2=10.0
       end
    3: begin
        coll1=20.0
        coll2=20.0
       end
    4: begin
        coll1=20.0
        coll2=40.0
       end
    5: begin
        coll1=40.0
        coll2=20.0
       end
    6: begin
        coll1=40.0
        coll2=40.0
       end
    7: begin
        coll1=40.0
        coll2=60.0
       end
    8: begin
        coll1=60.0
        coll2=40.0
       end
    9: begin
        coll1=60.0
        coll2=60.0
       end
    endcase
;----------------------------------- got the initial params

col1=coll1/60.0/57.269
col2=coll2/60.0/57.269
rmos=30./60./57.269
tmpd=col1^2+col2^2+4.*rmos^2
tmpn=(col1*col2)^2+(col1*rmos)^2+(col2*rmos)^2
dtheta=sqrt(tmpn/tmpd)

resol=fltarr(381,2)

for i=0,380 do begin
    en=i+1
    resol(i,0)=en          ;the energy array
    wavel=sqrt(81.805/en)
    stheta=wavel/(2.*dmono)
    if (abs(stheta) lt 1) then begin
       cot_theta=sqrt(1-stheta*stheta)/stheta
    endif else  begin
       cot_theta=0.0
    endelse
    de=2.*en*cot_theta*dtheta
    resol(i,1)=sqrt(de*de+1.2*1.2)
endfor

(*pstate).resfn=resol*8.06544      ; in wavenumbers!
widget_control,event.top,set_uvalue=pState,/no_copy
gaussfinalspectrum,event
end

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gaussfunGenWinDraw,event
; Event handler for plotting window.  Here we allow user to
; click with the left mouse button and
widget_control,event.top,get_uvalue = pState,/no_copy
if (n_elements(*(*pstate).freqptr) eq 0) then begin
  widget_control,event.top,set_uvalue = pState,/no_copy
  return
endif
xminrange=(*pstate).xminrange
xmaxrange=(*pstate).xmaxrange
yminrange=(*pstate).yminrange
ymaxrange=(*pstate).ymaxrange

case event.type of
0:  begin       ; button press
      (*pState).mouse = event.press
      if (*pState).mouse eq 1 then begin     ; what to do with a left click?
        widget_control,(*pState).plotWin1,draw_motion_events = 1
        (*pState).xbox[0] = event.x
        (*pState).ybox[0] = event.y
        wset,(*pState).winVis
        device,copy = [0,0,!d.x_size,!d.y_size,0,0,(*pState).winPix]
        empty
      endif
     if (*pState).mouse eq 4 then begin    ; what to do with a right click?
          (*pState).xminrange=min((*pState).freqspec)-10.0
          (*pState).xmaxrange=max((*pState).freqspec)+10.0
          (*pState).yminrange=min((*pState).fundamental)-(2.0)
      (*pState).ymaxrange=max((*pState).fundamental)+(2.0)
          if ((*pstate).plotover eq 1) then begin
          (*pState).yminrange=min((*pState).overtone)-(2.0)
          (*pState).ymaxrange=max(((*pState).overtone)+ (*pState).fundamental)+(2.0)
          endif
      widget_control,event.top,set_uvalue = pState,/no_copy
      gaussPlotcalc,event
    endif
end
1:  begin       ; button release

      if (*pState).mouse eq 1 then begin   ; what to do with a left mouse click?
        ; get the lower left hand corner of the zoom box in device coordinates
        xll = (*pState).xbox[0] < (*pState).xbox[1]
        yll = (*pState).ybox[0] < (*pState).ybox[1]
       ; calculate the width and height of the zoom box in device coordinates
        w = abs((*pState).xbox[1]-(*pState).xbox[0])
        h = abs((*pState).ybox[1]-(*pState).ybox[0])
        ; calculate the upper right corner of the zoom box
        xur = xll + w
        yur = yll + h
        ; convert the lower left coordinates from device to data
        ll = convert_coord(xll,yll,/device,/to_data)
        ; convert the upper right coordinates from device to data
        ur = convert_coord(xur,yur,/device,/to_data)

        (*pState).xminrange = ll[0]
        (*pState).xmaxrange = ur[0]
        (*pState).yminrange = ll[1]
        (*pState).ymaxrange = ur[1]

        wset,(*pState).winPix
        widget_control,event.top,set_uvalue = pState,/no_copy
       gaussPlotcalc,event
       widget_control,event.top,get_uvalue = pState,/no_copy
        wset,(*pState).winVis
        device,copy = [0,0,!d.x_size,!d.y_size,0,0,(*pState).winPix]
        empty

        widget_control,(*pState).plotWin1,draw_motion_events = 0
        (*pState).mouse = 0B
      endif
    end

2:  begin       ; mouse motion
    if (*pState).mouse eq 1 then begin ; what to do with a left mouse click?
      (*pState).xbox[1] = event.x  ; current mouse x-position
      (*pState).ybox[1] = event.y  ; current mouse y-position
      ; what are the box limits in device coordinates?
      xc = [(*pState).xbox[0],event.x,event.x,(*pState).xbox[0],(*pState).xbox[0]]
      yc = [(*pState).ybox[0],(*pState).ybox[0],event.y,event.y,(*pState).ybox[0]]
      wset,(*pState).winVis
      device,copy = [0,0,!d.x_size,!d.y_size,0,0,(*pState).winPix]
     plots,xc,yc,/device
      empty
    endif

    end
else:
endcase
widget_control,event.top,set_uvalue = pState,/no_copy
return


end

;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
pro gaussian98_event,event
    widget_control,event.top,get_uvalue=pState
    if (ptr_valid((pstate)) eq 0) then return
    WIDGET_CONTROL,EVENT.id,get_uvalue=uval
    freqchange=0


if (event.id eq ((*pState).scaleint)) then begin
       widget_control,(*pState).scaleint,get_value=newintscale
       scalefactor=(newintscale[0])
       (*pState).scalefactor=scalefactor
       widget_control,event.top,set_uvalue=pState
       if (n_elements(*(*pstate).iriptr) ne 0) then begin
       gaussfinalspectrum,event
       endif
       return

endif

if (event.id eq ((*pState).cwoversf)) then begin
       widget_control,(*pState).cwoversf,get_value=newscale
       oversf=(newscale[0])
       (*pState).oversf=oversf
       widget_control,event.top,set_uvalue=pState
       if (n_elements(*(*pstate).iriptr) ne 0) then begin
       gaussfinalspectrum,event
       endif
       return
endif

if (event.id eq ((*pState).bglinear)) then begin
       widget_control,(*pState).bglinear,get_value=newbg
       flatbg=(newbg[0])
       (*pState).flatbg=flatbg
       widget_control,event.top,set_uvalue=pState
       if (n_elements(*(*pstate).iriptr) ne 0) then begin
       gaussfinalspectrum,event
       endif
       return
endif

if (event.id eq ((*pState).dwfl)) then begin
       widget_control,(*pState).dwfl,get_value=newdw
       dwf=(newdw[0])
       (*pState).dwf=dwf
       widget_control,event.top,set_uvalue=pState
       if (n_elements(*(*pstate).iriptr) ne 0) then begin
       gaussfinalspectrum,event
       endif
       return
endif



if (event.id eq ((*pState).bgparam)) then begin
       widget_control,(*pState).bgparam,get_value=newbg
       linearbg=(newbg[0])
       (*pState).linearbg=linearbg
       widget_control,event.top,set_uvalue=pState
       if (n_elements(*(*pstate).iriptr) ne 0) then begin
       gaussfinalspectrum,event
       endif
       return
endif


if (event.id eq ((*pState).scalefreq)) then begin
       widget_control,(*pState).scalefreq,get_value=newfreqscale
         if (n_elements(*(*pstate).freqptr) eq 0) then return
       freq=(*(*pstate).freqptr)
       freqfactor=(newfreqscale[0])
       (*pState).freqfactor=freqfactor
       widget_control,event.top,set_uvalue=pState
       if (n_elements(*(*pstate).iriptr) ne 0) then begin
       gaussfinalspectrum,event
       endif
       return

endif

    case uval of
    'edatoms': begin
              if (n_elements((*((*pState).atomorig).listptr)) eq 0) then begin
                        print,'You must have atoms to edit them!'
                        return
              endif
              list=(*((*pState).atomorig).listptr)
              atom=(*((*pState).atomnew).newlistptr)
              ;atnum=(*((*pState).atomorig).numberptr)
              ;newatnum=(*((*pState).atomnew).newnumberptr)
              ixsection=(*((*pState).atomorig).ixsectionptr)
              IXS=(*((*pState).atomnew).newixsectionptr)
              mass=(*((*pState).atomnew).newmassptr)
              origmass=(*((*pState).atomorig).massptr)

              inatom={list:list,$
              origmass:origmass,$
              ixsection:ixsection}

              outatom={atom:atom,$
              mass:mass,$
              IXS:IXS}
              changed=0

              atli=gausseditatoms(group_leader=event.top,inatom,outatom,changed)

              (*((*pState).atomnew).newlistptr)=atli.atom
              (*((*pState).atomnew).newmassptr)=atli.mass
              (*((*pState).atomnew).newixsectionptr)=atli.IXS
              widget_control,event.top,set_uvalue=pState
              if (changed)  then begin
                 gaussEIGENVECTOR,event
              endif
              return
          end
'tomev': begin
          widget_control,event.id,get_value=val
          if (n_elements(*(*pstate).iriptr) eq 0) then begin
              return
          endif else if (val eq 'Freq. in wavenumber') then begin
              widget_control,event.id,set_value='Frequency is in meV'
              widget_control,event.top,/update
              (*pstate).freqspec=(*pstate).freqspec/8.06544
              (*(*pstate).freqptr)=(*(*pstate).freqptr)/8.06544
              (*pstate).resfn=(*pstate).resfn/8.06544
              freqchange=1
          endif else begin
              widget_control,event.id,set_value='Freq. in wavenumber'
              (*pstate).freqspec=(*pstate).freqspec*8.06544
              (*(*pstate).freqptr)=(*(*pstate).freqptr)*8.06544
              (*pstate).resfn=(*pstate).resfn*8.06544
              freqchange=1
          endelse
          widget_control,event.top,set_uvalue=pstate
         end
    'options': begin
          widget_control,event.id,get_value=val
          if (val eq 'Plot Options -->') then begin
              toolsbase=widget_base(event.top,/col)
              (*pstate).toolsbase=toolsbase
                                        a1toolsbase=widget_base(toolsbase,/row,frame=2)
              tool1=widget_button(a1toolsbase,value='Autoscale',event_pro='gaussautoscaleplot')
               tool2=cw_field(a1toolsbase,title = 'x min: ',/float,value = xminrange,$
                   uvalue = 'xminrange',/return_events,xsize = 10)
              tool3=cw_field(a1toolsbase,title = 'x max: ',/float,value = xmaxrange,$
                    uvalue = 'xmaxrange',/return_events,xsize = 10)
                tool4=cw_field(a1toolsbase,title = 'y min: ',/float,value = yminrange,$
                    uvalue = 'yminrange',/return_events,xsize = 10)
              tool5=cw_field(a1toolsbase,title = 'y max: ',/float,value = ymaxrange,$
                    uvalue = 'ymaxrange',/return_events,xsize = 10)
                   (*pstate).tool2=tool2
                   (*pstate).tool3=tool3
                   (*pstate).tool4=tool4
                  (*pstate).tool5=tool5




              widget_control,event.id,set_value='Less Options<--'
              widget_control,event.top,/update
          endif else begin
              widget_control,(*pstate).toolsbase,/destroy
              widget_control,event.id,set_value='Plot Options -->'
          endelse
            end
    'xminrange' : begin
              widget_control,event.id,get_value=xminrange
              (*pstate).xminrange=xminrange
              ;widget_control,event.top,set_uvalue=pState
                end
    'xmaxrange' : begin
              widget_control,event.id,get_value=xmaxrange
              (*pstate).xmaxrange=xmaxrange
          ;   widget_control,event.top,set_uvalue=pState
            end
    'yminrange' : begin
              widget_control,event.id,get_value=yminrange
              (*pstate).yminrange=yminrange
          ;   widget_control,event.top,set_uvalue=pState
            end
    'ymaxrange' : begin
              widget_control,event.id,get_value=ymaxrange
              (*pstate).ymaxrange=ymaxrange
              ;widget_control,event.top,set_uvalue=pState
            end

    'fund'  : begin
          if (event.select) then begin
          (*pState).plotfund = 1
          endif else begin
          (*pState).plotfund = 0
          endelse
          gaussautoscaleplot,event
          end
    'overtone'  : begin
          if (event.select) then begin
          (*pState).plotover = 1
          widget_control,event.top,set_uvalue=pState
          gaussfinalspectrum,event
          endif else begin
          (*pState).plotover = 0
          endelse
          gaussautoscaleplot,event
          end
    'data'  : begin
          if (event.select) then begin

          (*pState).plotdata = 1
          endif else begin
          (*pState).plotdata = 0
          endelse
          gaussautoscaleplot,event
          end
       else:
    endcase

widget_control,event.top,set_uvalue=pState

    if (n_elements(*(*pstate).spectrumptr) eq 0) then begin
       return
    endif else if (freqchange eq 1) then begin
       gaussautoscaleplot,event
    endif else  begin

       gaussplotcalc,event

    endelse
end

;;;;;;;;;;;;;;;;;;;;;;;;;;;;

pro gaussian98,group_leader = group_leader,datadirectory=datadirectory
;
; NB 11.6045 K =1meV and  1meV = 8.0655 cm-1, and 1 meV=0.2417 THz
;
; Extracts job info, and pertinant results (for me anyway).
;
;
;

if n_elements(group_leader) eq 0 then begin
    group_leader = 0L
endif
if n_elements(dataDirectory) eq 0 then dataDirectory = ''



  ; Get the current colors
  tvlct,r,g,b,/get
  device,get_decomposed = dc
  device,decomposed = 0
;loadct,0,/silent
;cd,current = thisDir
;cd,thisDir
!except = 0
;if n_elements(thisfile_path) eq 0 then thisfile_path = '/home/craig/idl/new/'


file_path=datadirectory


;isotopeinfo= ['H,1.01,80.27',$      ; name, mass, ixs
;'C,12.01,0.001',$
;'N,14.003,0.5',$
;'Ca,40.0,0.0']


; name, weigh AMU, COH Scat. Cross, Inc Scatt Cross, Total SCatt Cross, Abs Cross
isotopeinfo= [$
'H,1.00797,1.7568,80.26,82.02,0.3326',$
'D,2.014102,5.592,2.05,7.64,5.19E-4',$
'He,4.0026,1.34,0.0,1.34,0.00747',$
'Li,6.942,0.454,0.92,1.37,70.5',$
'Be,9.0122,7.63,0.0018,7.63,0.0076',$
'B,10.811,3.54,1.7,5.24,767.0',$
'C,12.01115,5.551,0.0010,5.551,0.0035',$
'N,14.0067,11.01,0.5,11.51,1.9',$
'O,15.9994,4.232,8.0E-4,4.232,1.9E-4',$
'F,18.9984,4.017,8.0E-4,4.018,0.0096',$
'Ne,20.183,2.62,0.0080,2.628,0.039',$
'Na,22.9898,1.66,1.62,3.28,0.53',$
'Mg,24.312,3.631,0.08,3.71,0.063',$
'Al,26.9815,1.495,0.0082,1.503,0.231',$
'Si,28.086,2.163,0.0040,2.167,0.171',$
'P,30.9738,3.307,0.0050,3.312,0.172',$
'S,32.064,1.0186,0.0070,1.026,0.53',$
'Cl,35.453,11.5257,5.3,16.8,33.5',$
'Ar,39.948,0.458,0.225,0.683,0.675',$
'K,39.102,1.69,0.27,1.96,2.1',$
'Ca,40.08,2.78,0.05,2.83,0.43',$
'Sc,44.956,19.0,4.5,23.5,27.5',$
'Ti,47.9,1.485,2.87,4.35,6.09',$
'V,50.942,0.0184,5.08,5.1,5.08',$
'Cr,51.996,1.66,1.83,3.49,3.05',$
'Mn,54.938,1.75,0.4,2.15,13.3',$
'Fe,55.847,11.22,0.4,11.62,2.56',$
'Co,58.9332,0.779,4.8,5.6,37.18',$
'Ni,58.71,13.3,5.2,18.5,4.49',$
'Cu,63.54,7.485,0.55,8.03,3.78',$
'Zn,65.37,4.054,0.077,4.131,1.11',$
'Ga,69.72,6.675,0.16,6.83,2.75',$
'Ge,72.59,8.42,0.18,8.6,2.2',$
'As,74.9216,5.44,0.06,5.5,4.5',$
'Se,78.96,7.98,0.32,8.3,11.7',$
'Br,79.909,5.8,0.1,5.9,6.9',$
'Kr,83.8,7.67,0.01,7.68,25.0',$
'Rb,85.47,6.32,0.5,6.8,0.38',$
'Sr,87.62,6.19,0.06,6.25,1.28',$
'Y,88.906,7.55,0.15,7.7,1.28',$
'Zr,91.22,6.44,0.02,6.46,0.185',$
'Nb,92.906,6.253,0.0024,6.255,1.15',$
'Mo,95.94,5.67,0.04,5.71,2.48',$
'Tc,99.0,5.8,0.5,6.3,20.0',$
'Ru,101.07,6.21,0.4,6.6,2.56',$
'Rh,102.905,4.34,0.3,4.6,144.8',$
'Pd,106.4,4.39,0.093,4.48,6.9',$
'Ag,107.87,4.407,0.58,4.99,63.3',$
'Cd,112.4,3.04,3.46,6.5,2520.0',$
'In,114.82,2.08,0.54,2.62,193.8',$
'Sn,118.69,4.871,0.022,4.892,0.626',$
'Sb,121.75,3.9,0.0070,3.9,4.91',$
'Te,127.6,4.23,0.09,4.32,4.7',$
'I,126.9044,3.5,0.31,3.81,6.15',$
'Xe,131.3,2.96,0.0,0,23.9',$
'Cs,132.905,3.69,0.21,3.9,29.0',$
'Ba,137.34,3.23,0.15,3.38,1.1',$
'La,138.91,8.53,1.13,9.66,8.97',$
'Ce,140.12,2.94,0.0010,2.94,0.63',$
'Pr,140.907,2.64,0.015,2.66,11.5',$
'Nd,144.24,7.43,9.2,16.6,50.5',$
'Pm,145.0,20.0,1.3,21.3,168.4',$
'Sm,150.35,0.422,39.0,39.0,5922.0',$
'Eu,151.96,6.57,2.5,9.2,4530.0',$
'Gd,157.25,29.3,151.0,180.0,49700.0',$
'Tb,158.925,6.84,0.0040,6.84,23.4',$
'Dy,162.5,35.9,54.4,90.3,994.0',$
'Ho,164.93,8.06,0.36,8.42,64.7',$
'Er,167.26,7.63,1.1,8.7,159.0',$
'Tm,168.934,6.28,0.1,6.38,100.0',$
'Yb,173.04,19.42,4.0,23.4,34.8',$
'Lu,174.97,6.53,0.7,7.2,74.0',$
'Hf,178.49,7.6,2.6,10.2,104.1',$
'Ta,180.948,6.0,0.01,6.01,20.6',$
'W,183.85,2.97,1.63,4.6,18.3',$
'Re,186.2,10.6,0.9,11.5,89.7',$
'Os,190.2,14.4,0.3,14.7,16.0',$
'Ir,192.2,14.1,0.0,14.0,425.0',$
'Pt,195.09,11.58,0.13,11.71,10.3',$
'Au,196.967,7.32,0.43,7.75,98.65',$
'Hg,200.59,20.24,6.6,26.8,372.3',$
'Tl,204.37,9.678,0.21,9.89,3.43',$
'Pb,207.19,11.115,0.0030,11.118,0.171',$
'Bi,208.98,9.148,0.0084,9.156,0.0338',$
'Po,210.0,0.0,0,00',$
'At,214.9986412,0.0,0,0,0',$
'Rn,222.0175705,0.0,0,12.6,0',$
'Fr,224.0,0.0,0,0,0',$
'Ra,226.0254,13.0,0.0,13.0,12.8',$
'Ac,227.0,0.0,0,0,0',$
'Th,2332.038,13.36,0.0,13.36,7.37',$
'Pa,232.5,10.4,0.1,10.5,200.6',$
'U,238.03,8.903,0.0050,8.908,7.57',$
'Np,237.0,14.0,0.5,14.5,175.9',$
'Pu,239.0522,0,0,0,0',$
'Am,242.0,8.7,0.3,9.0,75.3',$
'Cm,247.5,0.0,0,0,0']


isotopeinfoptr=ptr_new(isotopeinfo)

;some dummy file names
gfile=''
dfile=''
rotconstants=fltarr(3)   ; hold the rotational temperatures x,y,z in Kelvin
lowfrequencies=fltarr(9)    ; the low frequency vibrations  (floats)
job=strarr(20)       ; job type array (strings)
dof=0          ; degrees of freedom  (int)
maxforce=fltarr(2)   ; max force obtained, and threshold.
maxforceflag=''        ; did it meet the criteria
rmsforce=fltarr(2)   ; RMS force obtained, and threshold.
rmsforceflag=''        ; did it meet the criteria
maxdisp=fltarr(2)     ; max disp obtained, and threshold.
maxdispflag=''       ; did it meet the criteria
rmsdisp=fltarr(2)     ; RMS disp obtained, and threshold.
rmsdispflag=''       ; did it meet the criteria
numberatoms=0      ; how many atoms in the calculation


listptr=ptr_new(/allocate_heap)     ; atom list array (strings)
newlistptr=ptr_new(/allocate_heap)
numberptr=ptr_new(/allocate_heap)
newnumberptr=ptr_new(/allocate_heap)
ixsectionptr=ptr_new(/allocate_heap)
newixsectionptr=ptr_new(/allocate_heap)
massptr=ptr_new(/allocate_heap)
newmassptr=ptr_new(/allocate_heap)
initialorientationptr=ptr_new(/allocate_heap)

atomorig={listptr:listptr,$
         numberptr:numberptr,$
         ixsectionptr:ixsectionptr,$
         massptr:massptr}
atomnew={newlistptr:newlistptr,$
         newnumberptr:newnumberptr,$
         newixsectionptr:newixsectionptr,$
         newmassptr:newmassptr}

symmodeptr=ptr_new(/allocate_heap)
origfreqptr=ptr_new(/allocate_heap)
freqptr=ptr_new(/allocate_heap)
redmassptr=ptr_new(/allocate_heap)
frccptr=ptr_new(/allocate_heap)
iriptr=ptr_new(/allocate_heap)
ramanptr=ptr_new(/allocate_heap)
depolarptr=ptr_new(/allocate_heap)
vectorptr=ptr_new(/allocate_heap)
eigenyptr=ptr_new(/allocate_heap)   ; the eigenvector for each mode
spectrumptr=ptr_new(/allocate_heap) ; the final plot
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;   Want a big structure/pointer system to store the read file info...
;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
          readdata ={DatePtr:ptr_new(/allocate_heap),$
          ScanPtr:ptr_new(/allocate_heap),$
          MonitorPtr:ptr_new(/allocate_heap),$
          PrefactorPtr:ptr_new(/allocate_heap),$
          BasePtr:ptr_new(/allocate_heap),$
          TypePtr:ptr_new(/allocate_heap),$
          Collimationptr:ptr_new(/allocate_heap),$
          MosaicPtr:ptr_new(/allocate_heap),$
          XtalOrientPtr:ptr_new(/allocate_heap),$
          LatticeAPtr:ptr_new(/allocate_heap),$
          LatticeBPtr:ptr_new(/allocate_heap),$
          LatticeCPtr:ptr_new(/allocate_heap),$
          LatticeAlPtr:ptr_new(/allocate_heap),$
          LatticeBePtr:ptr_new(/allocate_heap),$
          LatticeGaPtr:ptr_new(/allocate_heap),$
          E_ZEROPtr:ptr_new(/allocate_heap),$
          DELTA_EPtr:ptr_new(/allocate_heap),$
          E_FINALPtr:ptr_new(/allocate_heap),$
          MON_D_SPACEPtr:ptr_new(/allocate_heap),$
          AN_D_SPACEPtr:ptr_new(/allocate_heap),$
          T_startPtr:ptr_new(/allocate_heap),$
          T_INCPtr:ptr_new(/allocate_heap),$
          HH_Ptr:ptr_new(/allocate_heap),$
          KK_Ptr:ptr_new(/allocate_heap),$
          LL_Ptr:ptr_new(/allocate_heap),$
          dH_Ptr:ptr_new(/allocate_heap),$
          dK_Ptr:ptr_new(/allocate_heap),$
          dL_Ptr:ptr_new(/allocate_heap),$
          HFIELD_Ptr:ptr_new(/allocate_heap),$
          QxPtr:ptr_new(/allocate_heap),$
          QyPtr:ptr_new(/allocate_heap),$
          QzPtr:ptr_new(/allocate_heap),$
          EnergyPtr:ptr_new(/allocate_heap),$
          TempPtr:ptr_new(/allocate_heap),$
          TimePtr:ptr_new(/allocate_heap),$
          TotalIntPtr:ptr_new(/allocate_heap)}
          readdataptr= ptr_new(readdata)


;;;;;;;;;;;;;;;;

scalefactor=1000.0            ;the 'calculated' intensity scalefactor
oversf=100.0                 ; scalefactor for overtones and comb.
Freqfactor=0.965              ;the 'calculated' freq scalefactor
flatbg=0.0                 ;background to be added
linearbg=0.0                 ;and a linear scaling with E
DWF=0.02                   ;debye-waller factor
d_cu=1.278                 ; Cu(220) d-spacing A
d_pg=3.354                 ; PG(002) d-spacing
resfn=fltarr(381,2)           ; resolution function to 300 meV
                        ;[0]=E and [1]=resolution
monochromator=['Cu','Pg','delta']
collimation=['10-10','10-20','20-10','20-20','20-40','40-20',$
                   '40-40','40-60','60-40','60-60']

xminrange=10.0
xmaxrange=250.0
yminrange=0.0
ymaxrange=1.0

fundamental=fltarr(2901)
overtone=fltarr(2901)
plotspec=fltarr(2901)
freqspec=findgen(2901)
freqspec=((freqspec/10)+10)*8.06544


tlb=widget_base(/col,title='FANS Gaussian',tlb_frame_attr=1,/base_align_left,mbar=bar)
filemenu = widget_button(bar,/menu,value='File')
void = widget_button(filemenu,value = 'Write Spectra as dave',event_pro='gaussdumpdave')
void = widget_button(filemenu,value = 'Write Spectra as ascii',event_pro='gaussdump')
void = widget_button(filemenu,value = 'Quit',event_pro='gaussquit')
filemenu2=widget_button(bar,value='Help')
void = widget_button(filemenu2,value = 'pdf help',event_pro='gausshelp')
void = widget_button(filemenu2,value = 'About',event_pro='gaussabout')

butbase1=widget_base(tlb,/row,/base_align_center,/frame)
choosefile=widget_button(butbase1,value='Choose Gaussian .log file',$
    event_pro='gaussChooseFile',xsize=160,tooltip='Load a LOG file')
choosedata=widget_button(butbase1,value='Choose expt. data file',$
    event_pro='gaussChooseData',xsize=160,tooltip='load any FANS file')
editatom=widget_button(butbase1,value='Edit atoms',$
    uval='edatoms',xsize=90,tooltip='change mass and ix of atoms')


monof=widget_droplist(butbase1,title="Monochromator",value=monochromator, $
       UValue=monochromator,Event_Pro='gaussmono_Events',xsize=160)
collimf=widget_droplist(butbase1,title="Collim",value=collimation, $
       UValue=collimation,Event_Pro='gausscoll_Events',xsize=120)
;writebutton=widget_button(butbase1,value='Write',$
;   event_pro='gaussdump',xsize=50,tooltip='Dump text files')

butbase3=widget_base(tlb,/row,/base_align_left)
butbase4=widget_base(butbase3,/col,/base_align_left,/frame)
;want to put a left pane of values and a right pane plot....

scaleint = cw_field(butbase4,title = 'Int. Scale factor: ',/float,value = scalefactor,$
                /return_events,xsize = 21)
scalefreq = cw_field(butbase4,title = 'Freq. Scale factor:',/float,value = freqfactor,$
                /return_events,xsize = 20)

butbase6=widget_base(butbase4,/row,/base_align_left)
bglinear = cw_field(butbase6,title = 'Bg, B0: ',value = flatbg,/float,$
                /return_events,xsize = 10)
bgparam = cw_field(butbase6,title = ',E*B1',value = linearbg,/float,$
                /return_events,xsize = 10)
dwfl = cw_field(butbase4,title = 'Debye-Waller factor:',value = dwf,/float,$
                /return_events,xsize = 18)
cwoversf = cw_field(butbase4,title = 'Comb. Scale fact:',value = oversf,/float,$
                /return_events,xsize = 20)

row1ne = widget_base(butbase4,/row, /NONEXCLUSIVE)
fundcheck = widget_button(row1ne, VALUE = 'Fundamentals', $
    UVALUE = 'fund',tooltip='plot the calculated fundamental')
overcheck = widget_button(row1ne, VALUE = 'Overtones', $
    UVALUE = 'overtone',tooltip='plot the calculated comb/overtones')
datacheck = widget_button(row1ne, VALUE = 'Data', $
    UVALUE = 'data',tooltip='plot the loaded experimental data')
widget_control, fundcheck, /SET_BUTTON


butbase7=widget_base(butbase4,/row,/base_align_center)
mevbutton=widget_button(butbase7,value='Freq. in wavenumber',uvalue='tomev',tooltip='Change energy scale')
morebutton=widget_button(butbase7,value='Plot Options -->',uvalue='options',tooltip='axis scaling')



winxsize=425 & winysize=278
butbase5=widget_base(butbase3,/row,/base_align_center)
plotwin1=widget_draw(butbase5,xsize=winxsize,ysize=winysize,$
       /button_events,event_pro = 'gaussfunGenWinDraw')


plotdata=0
plotover=0
plotfund=1



widget_control,tlb,/realize             ;put on screen

widget_control,plotwin1,GET_VALUE=winvis
window,/free,/pixmap,xsize = winxsize,ysize = winysize
winPix = !d.window

state = {file_path:file_path,$
        dump:0,$
    dumpdave:0,$
        mouse:0B,$
        gfile:gfile,$
    dfile:dfile,$
    readdataptr:readdataptr,$
    eigenyptr:eigenyptr,$
    spectrumptr:spectrumptr,$
    rotconstants:rotconstants,$
    lowfrequencies:lowfrequencies,$
    job:job,$
    dof:dof,$
    initialorientationptr:initialorientationptr,$
    isotopeinfoptr:isotopeinfoptr,$
    maxforce:maxforce,$
    maxforceflag:maxforceflag,$
    rmsforce:rmsforce,$
    rmsforceflag:rmsforceflag,$
    maxdisp:maxdisp,$
    maxdispflag:maxdispflag,$
    rmsdisp:rmsdisp,$
    rmsdispflag:rmsdispflag,$
    numberatoms:numberatoms,$
    atomorig:atomorig,$
    atomnew:atomnew,$
    symmodeptr:symmodeptr,$
    freqptr:freqptr,$
    origfreqptr:origfreqptr,$
    redmassptr:redmassptr,$
    frccptr:frccptr,$
    iriptr:iriptr,$
    ramanptr:ramanptr,$
    depolarptr:depolarptr,$
    vectorptr:vectorptr,$
    scalefactor:scalefactor,$
    freqfactor:freqfactor,$
    flatbg:flatbg,$
    linearbg:linearbg,$
    scaleint:scaleint,$
    scalefreq:scalefreq,$
    monochromator:monochromator,$
    collimation:collimation,$
    bglinear:bglinear,$
    bgparam:bgparam,$
    dwf:dwf,$
    dwfl:dwfl,$
    cwoversf:cwoversf,$
    oversf:oversf,$
    d_cu:d_cu,$
    d_pg:d_pg,$
    resfn:resfn,$
    monof:monof,$
    collimf:collimf,$
    mevbutton:mevbutton,$
    morebutton:morebutton,$
    xminrange:xminrange,$
    xmaxrange:xmaxrange,$
    yminrange:yminrange,$
    ymaxrange:ymaxrange,$
    plotdata:plotdata,$
    plotover:plotover,$
    plotfund:plotfund,$
    plotspec:plotspec,$
    fundamental:fundamental,$
    overtone:overtone,$
    freqspec:freqspec,$
    toolsbase:0L,$
    tool2:0L,$
    tool3:0L,$
    tool4:0L,$
    tool5:0L,$
    tlb:tlb,$
    plotwin1:plotwin1,$
    winvis:winvis,$
    winPix:winPix,$
        winxsize:winxsize,$
        winysize:winysize,$
    xbox:[0.0,1.0],$      ;device coordinates
    ybox:[0.0,1.0]};,$      ;device coordinates
;   inptr:inptr}

pState=ptr_new(state,/no_copy)

widget_control,tlb,set_uvalue=pState,/no_copy
xmanager,'gaussian98',tlb,cleanup = 'gausscleanup'

return
end
