; $Id: eis_moment__define.pro 405 2010-07-04 11:21:38Z viggoh $
;+
; NAME:
;       EIS_MOMENT__DEFINE
;
; PURPOSE:
;       EIS_MOMENT__DEFINE defines the class EIS_MOMENT. Objects of this
;       class contains intensity, velocity and line width data 
;       (0th, 1st and 2nd order intensity moments)
;       from the EIS instrument on SOLAR-B. The
;       EIS_MOMENT class inherits the superclass HW_DATA
;
; CATEGORY:
;       Hansteen/Wikstol Data analysis SW
;
; CALLING SEQUENCE:
;       The EIS_MOMENT__DEFINE procedure is not called directly. An
;       object of class EIS_INT_DATA is created with the following
;       statement: 
;                   eis_int_data = obj_new('eis_int_data', eis_data )
;
; INPUTS:
;       eis_data: Object of class EIS_DATA, containing the EIS data
;       from which the intensity data will be calculated and stored in
;       the instance data of the EIS_INT_DATA object.
;
; KEYWORD PARAMETERS:
;
;
; OUTPUTS:
;       Objects of type EIS_MOMENT 
;
; CALLS:
;
;
; COMMON BLOCKS:
;       
;
; PROCEDURE:
;       The procedure opens an object of class EIS_MOMENT. The
;       INIT-method of this object is called automatically. In the
;       INIT method, the method called CONSTRUCT is called. This
;       method uses the information from the EIS_DATA
;       object to perform moment calculations 
;       order to calculate intensity for all line window profiles. 
;       The moments are stored in the instance data called 'w' 
;       of EIS_MOMENT.
;  
; 
; RESTRICTIONS:
;
;
; MODIFICATION HISTORY:
;       04-Feb-2004: Viggo Hansteen. 
;       20-Feb-2004: Oivind Wikstol - Modified.
;       24-Feb-2004: Oivind Wikstol - Changed call procedure
;       14-Oct-2004: Oivind Wikstol - Changed data structure to
;                                     (nraster, yw, xw) where xw=3
;       15-Oct-2004: Oivind Wikstol - Added linelist: Pointer to index
;                                     list of processed lines. E.g. if total
;                                     number of lines are 25, but
;                                     moments have been calculated for
;                                     2 lines (3, 5), then *linelist
;                                     will be intarr(2), with
;                                     linelist[0]=3 and linelist[1]=5
;        3-Dec-2007: A. Gardini     - Added cleanup.
;        8-Feb-2008: A. Gardini     - Modified the save method. Moved
;                                     eis_fits_hdr_moment to construct.
;                                     Managed the missing of the error file.
;       16-Feb-2008: A. Gardini     - Gaussian fit computed by gauss_fit.
;       18-Mar-2008: A. Gardini     - Changed the column order of the 
;                                     datacube from [Y,X,m] to [m,Y,X]
;                                     for consistency with eis_data objects.
;                                     Added the mkfits method to save the
;                                     moments data in a level 2 file.
;       19-Mar-2008: A. Gardini     - Set the iwin keyword in save.
;        4-Jul-2008: A. Gardini     - Added HELP and INIT_EIS_MOMENTS_METHODS
;       11-Jul-2008: A. Gardini     - Set the group_leader variable.
;       19-Sep-2008: A. Gardini     - Introduced self.nfit .
;       12-Mar-2009: V. Hansteen    - Added aeffrelerr keyword, initiate 
;                                     calibration object as part of moment
;                                     object
;       30-Apr-2009: A. Gardini     - Added keywords to eis_cal
;        3-Aug-2009: A. Gardini     - Set self->setnwin to data->getnwin.
;-   
function eis_moment::init, data, iwin = iwin, fit = fit, $
                     group_leader = group_leader, $
                     lamtype = lamtype, line_def = line_def, $
                     aeffrelerr=aeffrelerr;, unit = unit
  self->init_eis_moment_methods ; initialization of the help method
; if n_elements(unit) eq 0 then unit = ['Counts', 'pixels', 'pixels!u2!n']
; self.unit = unit
  if n_elements(group_leader) gt 0 then self.group_leader = group_leader $
    else self.group_leader = -1
  if n_elements(lamtype) eq 0 then self.lam_type = 'wav' $
    else self.lam_type=lamtype
  if n_elements(line_def) ne 0 then self.line_def = 1 else self.line_def = 0 
  if n_elements(fit) eq 0 then self.fit='moment' else self.fit = fit
;  if n_elements(iwin) eq 0 then iwin = 0 else iwin = iwin
  self-> setcomment, 'EIS_moment'
  time_obs=(data->gethdr())->getdate_obs()
  slit_ind=data->getinfo('slit_ind')
  self.aux=ptr_new(obj_new('eis_aux'))
  self.cal=ptr_new(obj_new('eis_cal',time_obs = time_obs, slit_ind = slit_ind))
  if n_elements(aeffrelerr) ne 0 then begin
    *self.cal->setaeffrelerr,aeffrelerr
  endif
  if n_params() ge 1 then begin
    if n_elements(iwin) eq 0 then begin
      print, 'No data window specified. Processing all...'
      self-> construct, data, /all
    endif else begin
      print, 'Processing window ' +strtrim(iwin, 2) +'...'
      self-> construct, data, iwin = iwin
    endelse
  endif
  return, 1
end

pro eis_moment::construct, data, iwin = iwin, all = all, init=init
  if n_elements(init) eq 1 then return
  if n_params() lt 1 then begin
    print, 'eis_moment::construct, data, iwin = iwin, all = all'
    return
  endif
  if data->getfitslev() ne 1 then begin
    print, 'Data are not calibrated: moments cannot be computed.'
    return
  endif
  error=1
  if (data->geterr(0))[0] eq -1 then data->readerr,error=error
  if error eq -1 then begin
    print, 'Error data file is not available; moments cannot be computed.'
    return
  endif
  if not (data->getcalstat()).abs then $
    print,'WARNING. Absolute calibration was not performed. '+ $
          'The shown units are incorrect.'
  self.filename = data-> getfilename()
  self.datasource = data-> getdatasource()
  if ptr_valid(self.hdr) then ptr_free,self.hdr
  self.hdr = ptr_new(data-> gethdr())
 
  self.nraster = data->getnraster() ; & print,'nraster',self.nraster
  self.nexp=data->getnexp() ; & print,'nexp',self.nexp
  if ptr_valid(self.exp) then ptr_free,self.exp
  self.exp = ptr_new(fltarr(data->getnexp()))
  *self.exp = data->getexp()
  self.dx_size=data->getdx_size()
; self.dy_size=data->getdy_size()
  self.aux_data=data->getaux_data()
; self.nwin = data-> getnwin()
; nwin = self.nwin
  numwin = data-> getnwin()
  self-> setnwin, numwin
  if ptr_valid(self.linelist) then ptr_free,self.linelist
  self.linelist = ptr_new(intarr(numwin))
; if (self.fit eq 'moment') then self.nfit = 4 else self.nfit = 3
  self.nfit = 5 
  self.xs[0:numwin-1] = 0
  self.ys[0:numwin-1] = (data->getys())[0:numwin-1]
  self.xw[0:numwin-1] = self.nfit
  self.yw[0:numwin-1] = (data->getyw())[0]
  if n_elements(all) ne 0 then iwin=indgen(numwin)
  if n_elements(iwin) eq 0 then begin
    print, 'No line specified. Using default (0)'
    iwin = [0]
  endif
  iwin=data->getwindx(iwin)
  self.sit_and_stare=data->getsit_and_stare()

  nwin = n_elements(iwin)
  if ptr_valid(self.linelist) then ptr_free,self.linelist
  self.linelist = ptr_new(intarr(nwin))
  self.lam=ptr_new(fltarr(numwin,2))
  for iw = 0, nwin - 1 do begin
;    wd = data-> getvar(iwin[iw])
    (*self.linelist)[iw] = iwin[iw]
    (*self.lam)[iwin[iw],*]= $
       [min(data->getlam(iwin[iw])),max(data->getlam(iwin[iw]))]
    if ptr_valid(self.w[iwin[iw]]) then ptr_free,self.w[iwin[iw]]
    self.w[iwin[iw]] = $
;     ptr_new(fltarr(self.yw[iwin[iw]], self.nexp, self.xw[iwin[iw]]))
      ptr_new(fltarr(self.xw[iwin[iw]], self.yw[iwin[iw]], self.nexp)) 
    case self.fit of 
    'moment': begin
      self->moment,data,iwin[iw]
              end
    'gauss': begin
      self->gauss_fit,data,iwin[iw]
             end
    else: begin
          ok=dialog_message('Line fit option '+self.fit+' not implemented')
          end
    endcase
  endfor
; eis_fits_hdr_moment, self, *self.hdr, *self.linelist, numwin, data
end

pro eis_moment::cleanup
  ptr_free,self.linelist
  ptr_free,self.lam
  if ptr_valid(self.hdr) then ptr_free,self.hdr
  if ptr_valid(self.exp) then ptr_free,self.exp
  if ptr_valid(self.cal) then begin
    obj_destroy,*self.cal
    ptr_free,self.cal
  endif
  if ptr_valid(self.aux) then begin
    obj_destroy,*self.aux
    ptr_free,self.aux
  endif
  index = where(ptr_valid(self.w),ct)
  if ct gt 0 then ptr_free,self.w[index]
end

pro eis_moment::read, file, iwin = iwin, all = all, init=init
  if n_elements(init) eq 1 then return
  data = obj_new('eis_data', file)
  self-> construct, data, iwin = iwin, all = all
end

pro eis_moment::save, iwin=iwin, init=init
  if n_elements(init) eq 1 then return
  filename = self->getfilename()
  break_file,filename,disk,dir,name,ext
  name = str_replace(name,'_l1_','_l0_') ; if starting from a l1 file
  savefile = str_replace(name,'_l0_','_l2_')+'.fits'
  fxist =  file_test(savefile)
  if(fxist) then begin
    message =  'File already exists. Do you want to overwrite?'
    q =  dialog_message(message,/question,/center)
    if q eq 'No' then return
  endif
  if n_elements(iwin) eq 0 then iwin=indgen(self->getnwin())
  self->mkfits, filename, savefile, iwin
end

pro eis_moment::setlam_type, str, init=init
  if n_elements(init) eq 1 then return
  self.lam_type = strtrim(str, 2)
  return
end

function eis_moment::getlam_type,init=init
  if n_elements(init) eq 1 then return, -1
  return, self.lam_type
end

function eis_moment::getlinelist,init=init
  if n_elements(init) eq 1 then return, -1
  return, *self.linelist
end

pro eis_moment::setlam,iwin,lam
  (*self.lam)[iwin]=lam
end

function eis_moment::getlam,iwin
  return,(*self.lam)[iwin,*]
end

function eis_moment::mk_eis_map,win,mom=mom,int=int,vel=vel,wid=wid
  if keyword_set(wid) then mom=2
  if keyword_set(vel) then mom=1
  if keyword_set(int) then mom=0
  if n_elements(mom) eq 0 then mom=0
  iwin=(self->getwindx(win))[0]
  im=rotate(reform((self->getvar(iwin))[mom,*,*]),1)
  dx=(self->getxpos())[0]-(self->getxpos())[1]
  xc=(max(self->getxpos())+min(self->getxpos()))/2.
  yc=(max(self->getypos(iwin))+min(self->getypos(iwin)))/2.
  dy=((self->getypos())[1]-(self->getypos())[0])
  time=self->getinfo('DATE_OBS')
  eis_map=make_map(im,xc=xc,yc=yc,dx=dx,dy=dy,time=time)
  add_prop,eis_map,id='EIS '+self->getline_id(iwin),/replace
  dur=max(self->ti2tai(self->getti_2()))-min(self->ti2tai(self->getti_1()))
  add_prop,eis_map,dur=dur,/replace
  ang=pb0r(time,/arcsec)
  add_prop,eis_map,b0=ang[1]
  add_prop,eis_map,rsun=ang[2]
  return,eis_map
end


pro eis_moment::help
  help,self,/obj
  return
end

pro eis_moment::init_eis_moment_methods

  self->construct,/init
  self->read,/init
  self->save,/init
  self->setlam_type,/init
  a=self->getlam_type(/init)
  a=self->getlinelist(/init)
  self->corrvel,/init
  self->gauss_fit,/init
  self->mkfits,/init
  self->moment,/init
  self->tilt_correction,/init

  return
end

pro eis_moment__define
  struct={eis_moment, $
          p0:0.0, $
          fit:'gauss',$
          nfit:5,$
          lam_type:'pix', $
          line_def:0, $
          linelist:ptr_new(), $
          lam:ptr_new(), $
          group_leader:-1, $
          inherits eis_data}
end