pro calhrs_gimp,gimp,ihudl,udl,u1,u2,sdtype,ih,log
;+
;			hrs_gimp
;
; Routine to compute the geomagnetically induced image motion and correct
; the mapping function for the motion.
;
; CALLING SEQUENCE:
;	calhrs_gimp,gimp,ihudl,udl,u1,u2,ih,log
;
; INPUTS:
;	gimp - 2 element vector giving motion (diodes/gauss)
;		gimp(0) for detector 1
;		gimp(1) for detector 2
;	ihudl - udl headers (128 x nudl)
;	udl - udl array
;	u1 - vector of starting udl number for each science record
;	u2 - vector of ending udl number for each science record
;	sdtype - data type ('DDL' for direct downlink, 'ACC' for accum)
; INPUTS/OUTPUTS:
;	ih - header vectors for the science record.  Mapping function
;		info will be updated and motion value added.
;	log - processing log
;
; HISTORY:
;	version 1  D. Lindler  Feb. 1993
;-
;---------------------------------------------------------------------------

		VERSION = 1.0
;
; extract some useful info
;
	detector = ih(31)
	nudl = n_elements(ihudl)/128
	nsdp = n_elements(ih)/128
;
; get start and stop times for each readout -------------------------------
;
	if sdtype eq 'DDl' then begin
;
; direct down-link  (use science packet time as midpoint of data read)
;
		time = byte(ih(54:65),0,24,nsdp)
		time = string(time(0:23,*))
		juldates = jul_date(time)
	    end else begin
;
; accum mode (average first and last UDL times)
;
		pos = replicate(-1,nsdp)	;readout position
		juldates = dblarr(nsdp)		;juldates for each read
		nread = -1			;number of readouts
		u1last = -1			;begin udl for last bin
						;	processed		
		for i=0,nsdp-1 do begin
		    if (u1(i) ne -1) and (u2(i) ne -1) then begin ;udls present?
			if u1(i) ne u1last then begin          ;new readout?
				nread = nread+1
		    		jd1 = jul_date(string( $
					byte(ihudl(54:65,u1(i)),0,23)))	;start
		    		jd2 = jul_date(string( $
					byte(ihudl(54:65,u2(i)),0,23)))	;stop
				juldates(nread) = (jd1+jd2)/2.0
				u1last = u1(i)
			endif
			pos(i) = nread		;readout number for then bin
		    endif
		endfor
		juldates = juldates(0:nread)
	end		    
;
; compute epoch for earth's magnetic field
;
	date = date_conv(string(byte(ih(54:65),0,23)),'V')
	year = date(0) + date(1)/365.0		;epoch for mag. field
;
; conpute mag field in v1,v2,v3 coordinates
;
	st_geomag,log,year,juldates,bfield
;
; convert to hrs internal coordinates
;
	case detector of
		1: begin
			angle_a = 44.52		; angle about z axis
			angle_b = 0		; angle about x axis
			angle_c = 0		; angle about y axis
		   end
		2: begin
			angle_a = 45.0
			ANGLE_b = 0
			ANGLE_c = 0
		   end
	endcase
	rotate_3d,angle_a,3,bfield
	rotate_3d,angle_b,1,bfield
	rotate_3d,angle_c,2,bfield
;
; Rotate for the 17.6 degree EXB electron optical drift
;
	theta = 17.6/!radeg
	if n_elements(bfield) gt 3 then begin
		bx = -transpose(bfield(0,*))
		by = transpose(bfield(1,*))
	    end else begin
		bx = -bfield(0)
		by = bfield(1)
	end
	beff_x = cos(theta)*bx - sin(theta)*by
;
; compute gimp offset in sample (diode) units for each readout
;
	offsets = beff_x * gimp(detector-1)
;
; if ACCUM get offset for each bin
;
	if sdtype ne 'DDL' then begin
		offsets = offsets(pos>0)
		bad = where(pos lt 0,nbad)
		if nbad gt 0 then offsets(bad) = 0.0
	endif
;
; update sample positions in ih vectors and add gimp offset
;
	sample = float(ih(70:71,*),0,nsdp) - float(offsets)
	ih(70,0) = fix(sample,0,2,nsdp)
	ih(76,0) = fix(offsets,0,2,nsdp)
;
; update processing log
;
	if nbad gt 0 then begin
		hist = strarr(3)
		hist(2) = '    WARNING: unable to correct all data for GIMP;'+ $
			  ' Missing UDLs'
	end else hist = strarr(2)
	hist(0) = 'CALHRS_GIMP version '+string(version,'(f5.2)') + $
		': geomagnetic image motion correction'
	hist(1) = '    GIMP sensitivity = '+string(gimp(detector-1),'(F5.2)')+ $
		  ' Diodes/Gauss'
	if n_elements(log) gt 0 then sxaddhist,hist,log
	if !dump gt 0 then printf,!textunit,hist
return
end