PRO SIM_RASTER, UNIT, INPUT, HEADERS_ONLY=HEADERS_ONLY
;+
; Project     :	SOHO - CDS
;
; Name        :	SIM_RASTER
;
; Purpose     :	Simulates CDS telemetry for a raster
;
; Category    :	Telemetry, CDS
;
; Explanation :	This routine simulates the telemetry that one would expect from
;		a CDS raster, given the raster parameters.
;
; Syntax      :	SIM_RASTER, UNIT, INPUT
;
; Examples    :	
;
; Inputs      :	UNIT	= Logical unit number of file opened for write.  The
;			  telemetry packets are written into this file.
;		INPUT	= A structure variable containing the following
;			  parameters:
;
;				TIME     = Start time to use for first packet,
;					   in TAI format.  This value will be
;					   updated to reflect the times for the
;					   different packets.
;				COUNTER	 = Initial value of science data packet
;					   counter.  This number will be
;					   incremented for each packet written
;					   out.
;				SER_ID   = Series ID number
;				RAS_ID	 = Raster ID number.
;				RAS_VAR	 = Raster variation index.  The raster
;					   parameters are extracted from these
;					   two numbers.
;				SOLAR_X  = Solar X position
;				SOLAR_Y  = Solar Y position
;				STUDYCNT = Study counter
;				RAS_CNT  = Raster counter.  This is incremented
;					   by one each time this routine is
;					   called.
;
; Opt. Inputs :	None.
;
; Outputs     :	Certain parameters in the INPUT structure are updated by this
;		routine.
;
; Opt. Outputs:	None.
;
; Keywords    :	HEADERS_ONLY = If set, then only the raster and exposure
;			       headers are written out.  No dummy data is
;			       generated.
;
; Calls       :	
;
; Common      :	None.
;
; Restrictions:	The routine TM_DEFINE_STRUC must be called first to define all
;		the structures.
;
;		The output file must already be open for write access.
;
;		VDS readout orientation bit not yet accommodated.
;
; Side effects:	None.
;
; Prev. Hist. :	None.
;
; History     :	Version 1, 6-Jul-1995, William Thompson, GSFC
;		Version 2, 21-Aug-1995, William Thompson, GSFC
;			Added keyword HEADERS_ONLY
;
; Contact     :	WTHOMPSON
;-
;
	ON_ERROR, 2
;
;  Check the number of parameters.
;	
	IF N_PARAMS() NE 2 THEN MESSAGE, 'Syntax:  SIM_RASTER, UNIT, INPUT'
;
;  Get the raster definition.
;
	GET_RASTER, INPUT.RAS_ID, INPUT.RAS_VAR, RASTER
;
;  Not all compression schemes are supported.
;
	CASE RASTER.COMP_ID OF
		1: DUMMY = 0			;COPY is okay
		2: DUMMY = 0			;TRUNK is okay
		6: DUMMY = 0			;SUMLINE is okay
		7: DUMMY = 0			;SUMWIN is okay
		ELSE: RASTER.COMP_ID = 1	;Use COPY instead
	ENDCASE
;
;  Save the start time, and add the exposure time to the time.
;
	START_TIME = INPUT.TIME
	INPUT.TIME = INPUT.TIME + RASTER.EXPTIME
;
;  Form a generic packet header, and fill all the values in it.
;
	HEADER = {PACKETHDRDFN}
	HEADER.APID = '88A5'X
	HEADER.SOURCEPACKETCOUNTER = INPUT.COUNTER
	HEADER.DATAFIELDLENGTH = 293
	HEADER.TIME = TAI2OBT(INPUT.TIME)
	HEADER.STREAMPACKETCOUNTER = BYTE(INPUT.COUNTER)
	INPUT.COUNTER = INPUT.COUNTER + 1
;
;  Form the raster header packet.
;
	CASE RASTER.DETECTOR OF
		'G': BEGIN
			PACKET = {GISRASTERPKT}
			PACKET.HDR = HEADER
			PACKET.HDR.ID = '01'X
		     END
		'N': BEGIN
			PACKET = {VDSRASTERPKT}
			PACKET.HDR = HEADER
			PACKET.HDR.ID = '02'X
		     END
	ENDCASE
;
;  Next, fill all the parameters in the science part of the packet which is
;  common to both the GIS and NIS.
;
	PACKET.PKT.SCI.STUDYID = BYTE(INPUT.SER_ID)
	PACKET.PKT.SCI.COMPRESSION = RASTER.COMP_ID
	PACKET.PKT.SCI.DATAWINDOWID = RASTER.DW_ID
	PACKET.PKT.SCI.SEQUENCEID = INPUT.SER_ID
	PACKET.PKT.SCI.NUMMIRRORPOS = RASTER.NX
	PACKET.PKT.SCI.NUMSLITPOS = RASTER.NY
;
;  From the solar pointing, calculate the OPS L and R coordinates.
;
	GET_OPS_POS, INPUT.SOLAR_X, INPUT.SOLAR_Y, OPS_L, OPS_R
	PACKET.PKT.SCI.OPSLPOSITION = OPS_L
	PACKET.PKT.SCI.OPSRPOSITION = OPS_R
;
;  From the number of mechanism positions, and the step size, calculate the
;  initial mirror and slit positions.
;
	XSTEP = RASTER.XSTEP / 2
	YSTEP = RASTER.YSTEP
	NX = RASTER.NX > 1
	NY = RASTER.NY > 1
	MIRRORPOSITION = '80'X - FIX((NX - 1) / (2*XSTEP > 1))
	SLITPOSITION = -FIX((NY - 1) / (2*YSTEP > 1))
	PACKET.PKT.SCI.SLITPOSITION = SLITPOSITION
	PACKET.PKT.SCI.SLITNUMBER = RASTER.SLIT_NUM
	PACKET.PKT.SCI.MIRRORPOSITION = MIRRORPOSITION
;
;  The exposure time is complicated.  For NIS, it has to be an integer between
;  1 and 1023.  For GIS, it's always in units of 0.1 second.
;
	IF RASTER.DETECTOR EQ 'G' THEN BEGIN
	    EXPTIME = FIX(RASTER.EXPTIME*10)
	END ELSE BEGIN
	    EXPTIME = LONG(RASTER.EXPTIME*1000)
	    IF EXPTIME GT 1023 THEN EXPTIME = LONG(RASTER.EXPTIME*100)
	    IF EXPTIME GT 1023 THEN EXPTIME = LONG(RASTER.EXPTIME*10)
	    IF EXPTIME GT 1023 THEN EXPTIME = LONG(RASTER.EXPTIME)
	ENDELSE
	PACKET.PKT.SCI.EXPOSURETIME = EXPTIME
;
	PACKET.PKT.SCI.STUDYCOUNTER = INPUT.STUDYCNT
	PACKET.PKT.SCI.RASTERCOUNTER = INPUT.RAS_CNT
	PACKET.PKT.SCI.COMP_OPTION = RASTER.COMP_OPT
;
;  Calculate the window information.
;
	GET_DATAWIN, RASTER.DW_ID, DW
	N_WINDOWS = N_ELEMENTS(DW.WINS)
	PACKET.PKT.NWINDOWS = N_WINDOWS
	FOR I = 0,N_WINDOWS-1 DO BEGIN
	    CASE RASTER.DETECTOR OF
		'G': BEGIN
			PACKET.PKT.DATA_WIN(I).START = DW.WINS(I).WIN_DEF(0)
			PACKET.PKT.DATA_WIN(I).LENGTH = DW.WINS(I).WIN_DEF(2)
		     END
		'N': BEGIN
			PACKET.PKT.DATA_WIN(I).STARTX = DW.WINS(I).WIN_DEF(0)
			PACKET.PKT.DATA_WIN(I).STARTY = DW.WINS(I).WIN_DEF(1)
			PACKET.PKT.DATA_WIN(I).LENGTHX = DW.WINS(I).WIN_DEF(2)
			PACKET.PKT.DATA_WIN(I).LENGTHY = DW.WINS(I).WIN_DEF(3)
		     END
	    ENDCASE
	ENDFOR
;
;  Write out the raster header packet.
;
	SIM_WRT_PKT, UNIT, PACKET
;
;  Increment the packet time by 0.25 seconds to simulate the CDS data rate.
;
	INPUT.TIME = INPUT.TIME + 0.25
;
;  Step through the exposures.
;
	STARTMIRROR = MIRRORPOSITION
	EXPOSURENUMBER = 0
	FOR IY = 1,NY DO BEGIN
	    FOR IX = 1,NX DO BEGIN
;
;  Form the exposure header packet.
;
		PACKET = {EXPHEADERPKT}
		HEADER.SOURCEPACKETCOUNTER = INPUT.COUNTER
		HEADER.TIME = TAI2OBT(INPUT.TIME)
		HEADER.STREAMPACKETCOUNTER = BYTE(INPUT.COUNTER)
		INPUT.COUNTER = INPUT.COUNTER + 1
		PACKET.HDR = HEADER
		CASE RASTER.DETECTOR OF
			'G': PACKET.HDR.ID = '11'X
			'N': PACKET.HDR.ID = '12'X
		ENDCASE
;
;  From the OPS and mechanism positions, calculate the solar X and Y positions.
;
		OPS_POINT, FIX(OPS_L), FIX(OPS_R), SOLARX, SOLARY
		MECH_POINT, MIRRORPOSITION, SLITPOSITION, SOLARX, SOLARY
		PACKET.PKT.SOLAR_X = ROUND(SOLARX)
		PACKET.PKT.SOLAR_Y = ROUND(SOLARY)
;
;  Fill in the rest of the exposure header packet.
;
		PACKETCOUNTER = 0
		PACKET.PKT.PACKETCOUNTER = PACKETCOUNTER
		PACKET.PKT.STARTTIME = TAI2OBT(START_TIME)
		PACKET.PKT.EXPOSURENUMBER = EXPOSURENUMBER
		PACKET.PKT.OPSLPOSITION = OPS_L
		PACKET.PKT.OPSRPOSITION = OPS_R
		PACKET.PKT.SLITPOSITION = SLITPOSITION
		PACKET.PKT.MIRRORPOSITION = MIRRORPOSITION
		PACKET.PKT.STUDYCOUNTER = INPUT.STUDYCNT
		PACKET.PKT.RASTERCOUNTER = INPUT.RAS_CNT
;
;  Generate dummy data.
;
		IF NOT KEYWORD_SET(HEADERS_ONLY) THEN BEGIN
		    CASE RASTER.SLIT_NUM OF
			1:  BEGIN  &  SXWIDTH =  2  &  SYWIDTH =   2	& END
			2:  BEGIN  &  SXWIDTH =  4  &  SYWIDTH =   4	& END
			3:  BEGIN  &  SXWIDTH =  8  &  SYWIDTH =  50.8	& END
			4:  BEGIN  &  SXWIDTH =  2  &  SYWIDTH = 240	& END
			5:  BEGIN  &  SXWIDTH =  4  &  SYWIDTH = 240	& END
			6:  BEGIN  &  SXWIDTH = 90  &  SYWIDTH = 240	& END
			ELSE: BEGIN
				MESSAGE, 'Unrecognized slit number', /CONTINUE
				SXWIDTH = 2
				SYWIDTH = 2
				END
		    ENDCASE
		    IF RASTER.DETECTOR EQ 'G' THEN BEGIN
			GIS_DUMMY, ZONE=4, /NO_FILE, IMAGE=IMAGE, WLIST=1, $
				EXPOSE=RASTER.EXPTIME, SLIT=SXWIDTH,	$
				LENGTH=SYWIDTH
			FOR I = 0,N_WINDOWS-1 DO BEGIN
				I1 = DW.WINS(I).WIN_DEF(0)
				I2 = DW.WINS(I).WIN_DEF(2) + I1 - 1
				WIN = IMAGE(I1:I2)
				IF RASTER.COMP_ID GE 6 THEN WIN = TOTAL(WIN)
				IF I EQ 0 THEN DATA = WIN ELSE	$
					DATA = [DATA, WIN]
			ENDFOR
		    END ELSE BEGIN
			VDS_DUMMY, ZONE=4, /NO_FILE, IMAGE=IMAGE, WLIST=1, $
				EXPOSE=RASTER.EXPTIME, SLIT=SXWIDTH
			FOR I = 0,N_WINDOWS-1 DO BEGIN
				I1 = DW.WINS(I).WIN_DEF(0)
				I2 = DW.WINS(I).WIN_DEF(2) + I1 - 1
				J1 = DW.WINS(I).WIN_DEF(1)
				IF J1 GE 512 THEN J1 = 8 ELSE J1 = 128
				J2 = DW.WINS(I).WIN_DEF(3) + J1 - 1
				WIN = IMAGE(I1:I2,J1:J2)
				IF RASTER.COMP_ID EQ 6 THEN WIN = TOTAL(WIN,1)
				IF RASTER.COMP_ID EQ 7 THEN	$
					WIN = TOTAL(DOUBLE(WIN))
				IF I EQ 0 THEN DATA = WIN(*) ELSE	$
					DATA = [DATA, WIN(*)]
			ENDFOR
		    ENDELSE
;
;  Modify the data according to the compression scheme.
;
		    CASE RASTER.COMP_ID OF
			1:  DATA = FIX(DATA)
			2:  BEGIN
				TEMP = FIX(DATA) < 'FFF'X
				VDS_PACK, TEMP, DATA
			    END
			6:  DATA = LONG(DATA)
			7:  BEGIN
				TEMP = DATA
				MULT = 2.D0^32
				DATA = LONARR(2,N_ELEMENTS(TEMP))
				BIG = LONG(TEMP/MULT)
				TEMP = TEMP - BIG*MULT
				DATA(0,*) = LONG(TEMP)
				DATA(1,*) = BIG
				DATA = DATA(*)
			    END
		    ENDCASE
;
;  Calculate the total number of bytes to write out, and convert the data to a
;  byte array.
;
		    CASE DATATYPE(DATA,2) OF
			2:  N_BYTES = 2*N_ELEMENTS(DATA)
			3:  N_BYTES = 4*N_ELEMENTS(DATA)
		    ENDCASE
		    DATA = BYTE(DATA, 0, N_BYTES)
;
;  Store the first 252 bytes of the data in the exposure header.
;
		    TEMP = BYTARR(252)
		    TEMP(0) = DATA(0:(251<(N_BYTES-1)))
		    PACKET.PKT.SCIENCE = TEMP
		    N_BYTES = N_BYTES - 252
		    IF N_BYTES GT 0 THEN DATA = DATA(252:*) ELSE	$
			PACKET.PKT.PACKETCOUNTER = '8000'X  OR		$
				PACKET.PKT.PACKETCOUNTER
		END ELSE N_BYTES = 0
;
;  Write out the exposure header packet.
;
		SIM_WRT_PKT, UNIT, PACKET
;
;  Increment the packet time by 0.25.
;
		INPUT.TIME = INPUT.TIME + 0.25
;
;  Keep writing out science packets until there is no more data left.
;
		WHILE N_BYTES GT 0 DO BEGIN
			PACKET = {SCIPACKETPKT}
			HEADER.SOURCEPACKETCOUNTER = INPUT.COUNTER
			HEADER.TIME = TAI2OBT(INPUT.TIME)
			HEADER.STREAMPACKETCOUNTER = BYTE(INPUT.COUNTER)
			INPUT.COUNTER = INPUT.COUNTER + 1
			PACKET.HDR = HEADER
			CASE RASTER.DETECTOR OF
				'G': PACKET.HDR.ID = '21'X
				'N': PACKET.HDR.ID = '22'X
			ENDCASE
			PACKETCOUNTER = PACKETCOUNTER + 1
			PACKET.PKT.PACKETCOUNTER = PACKETCOUNTER
			PACKET.PKT.EXPOSURENUMBER = EXPOSURENUMBER
			PACKET.PKT.STUDYCOUNTER = INPUT.STUDYCNT
			PACKET.PKT.RASTERCOUNTER = INPUT.RAS_CNT
;
			TEMP = BYTARR(276)
			TEMP(0) = DATA(0:(275<(N_BYTES-1)))
			PACKET.PKT.SCIENCE = TEMP
			N_BYTES = N_BYTES - 276
			IF N_BYTES GT 0 THEN DATA = DATA(276:*)	ELSE	$
				PACKET.PKT.PACKETCOUNTER = '8000'X  OR	$
					PACKET.PKT.PACKETCOUNTER
;
;  Write out the science packet.
;
			SIM_WRT_PKT, UNIT, PACKET
;
;  Increment the packet time by 0.25.
;
			INPUT.TIME = INPUT.TIME + 0.25
		ENDWHILE			
;
;  Compare the next packet time to the start time plus the exposure time and
;  use whichever is larger.
;
		START_TIME = START_TIME + RASTER.EXPTIME
		IF INPUT.TIME GT START_TIME THEN START_TIME = INPUT.TIME
		IF START_TIME GT INPUT.TIME THEN INPUT.TIME = START_TIME
;
;  Increment the mirror and slit positions, and do the next exposure.
;
		EXPOSURENUMBER = EXPOSURENUMBER + 1
		MIRRORPOSITION = MIRRORPOSITION + XSTEP
	    ENDFOR
	    MIRRORPOSITION = STARTMIRROR
	    SLITPOSITION = SLITPOSITION + YSTEP
	ENDFOR
;
;  Increment the raster counter.
;
	INPUT.RAS_CNT = INPUT.RAS_CNT + 1
;
	RETURN
	END