function AMIE_COMBINE, INDIR = inDir, INFILE = inFile, ORIENTATION = orientation, $
                       AUTOSCALE = autoscale, silent = silent, $
                       integer = integer

;+
; $Log: AMIE_COMBINE.PRO,v $
; Revision 1.2  2008/02/14 12:06:39  bgrieger
; Changed default behaviour to return 4 byte floats and implemented option
; to return 2 byte integers.
; A few revisions had been lost when the repository was transfered from
; ESTEC to ESAC. Now all revisions should be back in.
; Modified Files:
; 	AMIE_COMBINE.PRO
;
; Revision 1.5  2007/07/04 10:15:53  bgrieger
;    V1.7 - 2007 Jul 04, 12h00m, bg  - /continue when a file is not found
;                                      (amie_combine does not stop but returns)
;                                      but do not perform /autoscale in this
;                                      case (avoids arithmetic overflow).
; Modified Files:
; 	AMIE_COMBINE.PRO
;
; Revision 1.4  2007/06/19 10:39:49  bgrieger
; Added passing through of the /SILENT option.
; Modified Files:
; 	AMIE_COMBINE.PRO READ_AMIE.PRO
;
; Revision 1.3  2007/06/13 17:10:16  bgrieger
; AMIE_COMBINE:
; ;   V1.5 - 2007 Jun 13, 18h00m, bg  - Do not assume anymore that the filter
; ;                                     number is character number 7 in the file
; ;                                     name, but use stregex to get its position
; ;                                   - Deleted a few lines which were commented
; ;                                     out anyway.
; ;
; READ_AMIE:
; ;    2007 Jun 13, V1.7, bg: Treat DATA_QUALITY_ID=-1 by first calling HEADPDS and not
; ;                       calling READPDS if this is the case. Otherwise non existing
; ;                       data would cause an error.
; ;
; Modified Files:
; 	AMIE_COMBINE.PRO READ_AMIE.PRO
;
; Revision 1.2  2007/05/28 10:04:04  bgrieger
; ;   V1.4 - 2007 May 21, 12h00m, bg  - Fixed errors in the AUTOSCALE option:
; ;                                       - corrected indices for computation of
; ;                                         mean and variance of FeH_Y filter
; ;                                       - introduced float array for the image
; ;                                         to deal with negative scaled values
; ;                                         that may occur under extreme
; ;                                         illumination conditions
; ;                                   - Deleted the writing of image data to a
; ;                                     file (which had been introduced for
; ;                                     testing)
; ;                                   - Deleted double declaration of image array
; ;                                   - Set initial image array really to grey
; Modified Files:
; 	AMIE_COMBINE.PRO
;
; Revision 1.1  2007/05/25 14:24:25  smartine
; 2007.05.25, SM: Added AMIE_COMBINE.PRO and READ_AMIE.PRO
;
; Revision 1.1  2007/05/25 06:46:14  smartine
; 2007.05.25, SM:  First version of this file
;
;
; NAME:
;   AMIE_COMBINE
;
; PURPOSE:
;   Reads eight AMIE label and image files containing the individual filter images
;   and generates one full-frame (1024 x 1024) image from it. This is returned.
;
; CATEGORY:
;   AMIE STANDARD software routines
;
; CALLING SEQUENCE:
;   image = AMIE_COMBINE ()
;
; INPUTS:
;   n/a
;
; OPTIONAL INPUT PARAMETERS:
;   INDIR:  Set this to a string variable denoting a valid path. If set, the program will
;           look in this path for input files.
;
;   INFILE: File name to be read.
;
;   ORIENTATION: Possible values are 'CCD_FRAME' or 'REAL_WORLD'. This keyword
;                specifies the orientation of the image. If this keyword is not
;                present, 'REAL_WORLD' will be used.
;                See also S1-AMI-RSSD-TN-001/1d, Figure 11.
;
;                CCD_FRAME will give the following orientation:
;
;                 y    ^
;                  CCD |
;                      |
;                      |
;                      |
;                      +-------> x
;                                 CCD
;
;                REAL_WORLD will give the following orientation:
;
;                      +-------> x
;                      |          CCD
;                      |
;                      |
;                      |
;                      v
;                       y
;                        CCD
;
;                This view corresponds to what an observer would see if he sat
;                on the spacecraft with his feet mounted on the side plane of
;                Smart-1 just as AMIE, looking in the direction of the line of
;                sight of the camera.
;
;   AUTOSCALE: If set to 1, all filter images are scaled to the same
;              mean and variance. As standard mean and variance either
;              the values of the CLEAR filter or of the FeH_Y filter are
;              used, depending on which of the two has the higher
;              variance.
;
;   INTEGER:   If set to one, the combined image is returned as 2 byte
;              integer (the original behaviour). As default, the image
;              is returned as 4 byte float.
;
; OUTPUTS:
;   The return value of the function is a full-frame (1024 x 1024) image
;
; OPTIONAL OUTPUT PARAMETERS:
;   None.
;
; COMMON BLOCKS:
;   None.
;
; SIDE EFFECTS:
;   Hopefully none.
;
; RESTRICTIONS:
;   This routine does not do consistency checking, i.e. if the file name of a filter
;   image is wrong the image will be messed up.
;
;   If one filter is missing, that part of the image will be black and a warning
;   will be issued in the message window of IDL. /AUTOSCALE will not work.
;
;   Only tested on Windows and PC-Linux systems, no guarantee for proper functioning
;   on other platforms.
;
; EXAMPLES:
;   The simplest way of calling AMIE_COMBINE is by typing
;
;   IDL> IMAGE = AMIE_COMBINE ()
;
;   This will open a DIALOG_PICKFILE () menu. Select an image file. The program will
;   search for the corresponding filter images taken at the same time and
;   the program will combine them to one full-frame image.
;
;
; PROCEDURE:
;
;   If INDIR is set, it will go to the directory and let the user pick a file. If
;   INFILE is set, it will just take that one.
;   According to the EAICD, p. 18, the name looks like this:
;   AMI_MPTF_XXXXXX_NNNNN_TTTTT.IMG, where 'F' is the filter number. All other
;   parts of the name must be the same if the images go together. So we do:
;
;   (a) generate a grey image, 1024 x 1024
;   (b) take the first filename
;   (c) find the corresponding ones
;   (d) put them on top of the grey image
;   (e) scale all filter images to similar brightness (if requested)
;   (f) return the image
;
;   The orientation of the image is as shown below (NOTE: The shown orientation
;   corresponds to the 'REAL_WORLD' view - for ORIENTATION = 'CCD_FRAME', the
;   image would be flipped horizontally, i.e. the +y_CCD axis would go to the
;   top).
;
;   The numbers denote the filter numbers. The readout register would be on top,
;   read out to the left. This orientation is such that the image is right side
;   up, i.e. if a human observer would look along the line-of-sight of the camera,
;   he/she would have the same view.
;
;   +------+------+------+------+----> x
;   |      |      |      |      |       CCD
;   |  (5) |      |      |      |
;   |      | (8)  |      |      |
;   +------+      | (6)  | (7)  |
;   |             |      |      |
;   |   (4)       |      |      |
;   |             |      |      |
;   +-------------+------+------+
;   |             |             |
;   |     (2)     |             |
;   |             |             |
;   +-------------+     (3)     |
;   |             |             |
;   |     (1)     |             |
;   |             |             |
;   +-------------+-------------+
;   |
;   v
;    y
;     CCD
;
; MODIFICATION HISTORY:
;   V0.1 - 2004 Sep 01, 13h00m, dvk: First go, derived from AMIE_EXTRACT
;   V0.2 - 2004 Sep 01, 19h00m, dvk: continue
;   V1.0 - 2004 Sep 01, 21h00m, dvk - NOTE: not properly tested yet!
;   V1.1 - 2005 Sep 27, 13h00m, dvk - updating comments, fixing image order
;   V1.2 - 2007 May 16, 11h20m, dvk, bg - in the text above, it said 'readout to
;                                         the right' - this is wrong, it has to say
;                                         'readout to the left'. Fixed.
;                                       - Delete the rotation of the individual images
;                                       - Add the keyword ORIENTATION
;                                       - Add tested on PC-Linux in
;                                         the text above
;   V1.3 - 2007 May 21, 16h30m, bg  - Add the AUTOSCALE option
;   V1.4 - 2007 May 21, 12h00m, bg  - Fixed errors in the AUTOSCALE option:
;                                       - corrected indices for computation of
;                                         mean and variance of FeH_Y filter
;                                       - introduced float array for the image
;                                         to deal with negative scaled values
;                                         that may occur under extreme
;                                         illumination conditions
;                                   - Deleted the writing of image data to a
;                                     file (which had been introduced for
;                                     testing)
;                                   - Deleted double declaration of image array
;                                   - Set initial image array really to grey
;   V1.5 - 2007 Jun 13, 18h00m, bg  - Do not assume anymore that the filter
;                                     number is character number 7 in the file
;                                     name, but use stregex to get its position
;                                   - Deleted a few lines which were commented
;                                     out anyway.
;   V1.6 - 2007 Jun 19, 13h00m, bg  - Added passing the /silent option to
;                                     READ_AMIE
;   V1.7 - 2007 Jul 04, 12h00m, bg  - /continue when a file is not found
;                                     (amie_combine does not stop but returns)
;                                     but do not perform /autoscale in this
;                                     case (avoids arithmetic
;                                     overflow).
;   V1.8 - 2008 Feb 14, 13h00m, bg  - Changed default behaviour to
;                                     return 4 byte float and
;                                     implemented option to return 2
;                                     byte integer.
;
; (c) 2004-2007, Detlef V. Koschny, ESA/RSSD
;-

;------------------------------------------------------------------------------
; Perform some checks on keywords (Label is handled later)
;------------------------------------------------------------------------------
if n_elements (INDIR) GT 0 then begin
    if file_test (INDIR, /dir) EQ 0 then $
       message, 'Input directory "' + INDIR + '" does not exist.'
    ;if not there, add the path separator at end ("\" for Windows)
    if strmid (InDir, strlen (InDir), 1) NE PATH_SEP () then $
                                                   InDir = InDir + PATH_SEP ()
endif

if n_elements (INFILE) GT 0 then begin
    if n_elements (INDIR) NE 1 then $
       message, 'If INFILE exists, INDIR must exist too'
    if file_test (INDIR + INFILE) EQ 0 then $
       message, 'File "' + INFILE + '" does not exist in dir ' + INDIR + '.', $
                /continue
endif

if keyword_set (ORIENTATION) then orientation = ORIENTATION else $
                                  orientation = 'REAL_WORLD'
if (orientation NE 'REAL_WORLD') and (orientation NE 'CCD_FRAME') then $
                       message, '[ERROR] Wrong value for keyword "ORIENTATION"'

if keyword_set (AUTOSCALE) then autoscale = autoscale else $
                                autoscale = 0
if (autoscale NE 0) and (autoscale NE 1) then $
                       message, '[ERROR] Wrong value for keyword "AUTOSCALE"'

if keyword_set (silent) then silent = silent else silent = 0

if keyword_set( integer ) then begin
    integer = integer
endif else begin
    integer = 0
endelse
if integer ne 0 then integer = 1

;------------------------------------------------------------------------------
; Defining arrays, constants, etc.
;------------------------------------------------------------------------------
filelist = strarr (1)
filter_number = strarr (9)
filter_name = strarr (9)
record_bytes = strarr (9)
file_records = strarr (9)
label_records = strarr (9)
center_wavelength = strarr (9)
bandwidth = strarr (9)
x_size = intarr (9)
y_size = intarr (9)
x_start = intarr (9)
y_start = intarr (9)
filename = strarr (1)

grey = 250
if integer eq 1 then begin
    image = uintarr (1024, 1024)
endif else begin
    image = fltarr (1024, 1024)
endelse
image_flt = fltarr (1024, 1024)

;------------------------------------------------------------------------------
; Defining the special things for the individual filter images
; Each filter gets a number, 0 is the full image. This follows the EAICD
; (S1-AMIE-RSSD-IF-001)
;------------------------------------------------------------------------------

;Image positions in full frame and size
;--------------------------------------
x_size (0) = 1024 & y_size (0) = 1024 & x_start (0) =    0 & y_start (0) =    0
x_size (1) = 512  & y_size (1) =  256 & x_start (1) =    0 & y_start (1) =    0
x_size (2) = 512  & y_size (2) =  256 & x_start (2) =    0 & y_start (2) =  256
x_size (3) = 512  & y_size (3) =  512 & x_start (3) =  512 & y_start (3) =    0
x_size (4) = 512  & y_size (4) =  256 & x_start (4) =    0 & y_start (4) =  512
x_size (5) = 256  & y_size (5) =  256 & x_start (5) =    0 & y_start (5) =  768
x_size (6) = 256  & y_size (6) =  512 & x_start (6) =  512 & y_start (6) =  512
x_size (7) = 256  & y_size (7) =  512 & x_start (7) =  768 & y_start (7) =  512
x_size (8) = 256  & y_size (8) =  512 & x_start (8) =  256 & y_start (8) =  512

;-------------------------------------------------------------------------------
; (a) generate a grey image, 1024 x 1024
;-------------------------------------------------------------------------------
image (*, *) = grey

;-------------------------------------------------------------------------------
; (b) get the name of one of the filter images, depending on keywords set
;-------------------------------------------------------------------------------
fileOk = 0
if n_elements (inDir) NE 1 then begin
  inFile = dialog_pickfile (get_path = inDir, filter = "*.img", /read)

  ;to be compatible with the rest of the program, we need to reduce inFile
  ;to the file name only.
  inFile = file_basename (inFile)
  fileOk = 1
endif

if (n_elements (inFile) NE 1) AND NOT (fileOk) then begin
  if n_elements (inDir) then begin
       inFile = dialog_pickfile (path = inDir, filter = "*.img", /read)
  endif else inFile = dialog_pickfile (filter = "*.img", /read)
  inFile = file_basename (inFile)   ;see above
endif

;-------------------------------------------------------------------------------
; (c) read all eight files  and (d) put on top of the grey image
;-------------------------------------------------------------------------------
;Get the position of the filter number in the filename
ipos = stregex( inFile, '^[A-Z_]+([1-8])_', /SUBEXPR )
if ipos[1] EQ -1 then begin
  message, 'ERROR: Filter number cannot be found in filename ' + inFile
endif
missing = 0
for i = 1, 8 do begin
  ;strcompress is needed, otherwise there are leading spaces when converting
  ;i to string.
  strput, inFile, strcompress (string (i), /remove_all), ipos[1]
  if file_test (inDir + inFile) EQ 0 then begin
    message, 'WARNING: File ' + inFile + ' is missing.', /continue
    missing = 1
  endif else begin
    ;file does exist, read in
    filterimage = READ_AMIE (inDir + inFile, label, $
        orientation = 'REAL_WORLD', silent = silent)
    ;put in right place - coordinates come from EAICD table
    image (x_start (i) : x_start (i) + x_size (i) - 1, $
           y_start (i) : y_start (i) + y_size (i) - 1) = filterimage.image1
  endelse
endfor

;-------------------------------------------------------------------------------
; (e) scale all filter images to the same mean and variance
;-------------------------------------------------------------------------------
if autoscale eq 1 and missing eq 0 then begin

   ; this strange back and forth to avoid type convertion of image_flt to
   ; integer:
   image_flt = -1000. * image
   image_flt = image_flt/(-1000.)

   ; mean and variance of the clear filter (avoiding bright left and
   ; upper edge):
   meanclear = mean( reform( image_flt[528:1023,0:495], 246016 ) )
   varclear = variance( reform( image_flt[528:1023,0:495], 246016 ) )

   ; mean and variance of the FeH_Y filter:
   meanfilter = mean     ( reform( image_flt[0:511,512:767], 131072 ) )
   varfilter  = variance ( reform( image_flt[0:511,512:767], 131072 ) )
   meanfehy = meanfilter
   varfehy  = varfilter

   ; make the one with the larger variance the standard for all others
   ; (stored in the -clear variables in any case) and scale the CLEAR filter
   ; if it is not the standard:
   if varfehy gt varclear then begin
      image_flt[512:1023,0:511] = meanfehy + $
         ( image_flt[512:1023,0:511] - meanclear ) * sqrt( varfehy / varclear )
      meanclear = meanfehy
      varclear  = varfehy
   endif

   ; upper and lower limit for the combined image:
   ulim = meanclear + 3.*sqrt(varclear)
   llim = meanclear - 3.*sqrt(varclear)
   if llim lt 0 then llim = 0

   ; scale FeH_X (avoiding bright upper edge) and FeH_Y
   meanfilter = (2./3.) * meanfilter + (1./3.) * $
      mean     ( reform( image_flt[256:511,768:991], 57344 ) )
   varfilter  = (2./3.) * varfilter  + (1./3.) * $
      variance ( reform( image_flt[256:511,768:991], 57344 ) )
   image_flt[0:511,512:767] = meanclear + $
      ( image_flt[0:511,512:767] - meanfilter ) * sqrt( varclear / varfilter )
   image_flt[256:511,768:1023] = meanclear + $
      ( image_flt[256:511,768:1023] - meanfilter ) * sqrt( varclear / varfilter )

   ; scale LASER (avoiding bright upper edge):
   meanfilter = mean     ( reform( image_flt[0:255,768:991], 57344 ) )
   varfilter  = variance ( reform( image_flt[0:255,768:991], 57344 ) )
   image_flt[0:255,768:1023] = meanclear + $
      ( image_flt[0:255,768:1023] - meanfilter ) * sqrt( varclear / varfilter )

   ; scale FeL_X (avoiding bright lower and upper edges):
   meanfilter = mean     ( reform( image_flt[512:767,528:991], 118784 ) )
   varfilter  = variance ( reform( image_flt[512:767,528:991], 118784 ) )
   image_flt[512:767,512:1023] = meanclear + $
      ( image_flt[512:767,512:1023] - meanfilter ) * sqrt( varclear / varfilter )

   ; scale VIS_X (avoiding bright lower and upper edges):
   meanfilter = mean     ( reform( image_flt[768:1023,528:991], 118784 ) )
   varfilter  = variance ( reform( image_flt[768:1023,528:991], 118784 ) )
   image_flt[768:1023,512:1023] = meanclear + $
      ( image_flt[768:1023,512:1023] - meanfilter ) * sqrt( varclear / varfilter )

   ; scale FeL_Y (avoiding bright right edge):
   meanfilter = mean     ( reform( image_flt[0:495,256:511], 126976 ) )
   varfilter  = variance ( reform( image_flt[0:495,256:511], 126976 ) )
   image_flt[0:511,256:511] = meanclear + $
      ( image_flt[0:511,256:511] - meanfilter ) * sqrt( varclear / varfilter )

   ; scale VIS_Y (avoiding bright right edge):
   meanfilter = mean     ( reform( image_flt[0:495,0:255], 126976 ) )
   varfilter  = variance ( reform( image_flt[0:495,0:255], 126976 ) )
   image_flt[0:511,0:255] = meanclear + $
      ( image_flt[0:511,0:255] - meanfilter ) * sqrt( varclear / varfilter )

   ; apply upper and lower limit:
   for i = 0, 1023 do begin
      for j = 0, 1023 do begin
         if ( image_flt[i,j] gt ulim ) then image_flt[i,j] = ulim
         if ( image_flt[i,j] lt llim ) then image_flt[i,j] = llim
      endfor
   endfor
   
   if integer eq 1 then begin
       image = round( image_flt )
   endif else begin
       image = image_flt
   endelse

endif

;-------------------------------------------------------------------------------
; (f) return the image
;-------------------------------------------------------------------------------
; If orientation EQ 'CCD_FRAME', flip horizontally
if orientation EQ 'CCD_FRAME' then image = rotate (rotate (image, 1), 4)

return, image

end ;function AMIE_COMBINE