Version 0.0

interpol_param.pro

+function interpol_param, hk_table, par_name, time
+	i = where(hk_table.par_name eq par_name, count)
+	hk_time = dblarr(count)
+	for j = 0, count - 1 do begin
+		times = strsplit(hk_table.gen_time[i[j]], ":", /extract)
+		hk_time[j] = double(times[0]) + double(times[1])/65536.0D
+	endfor 
+	value = interpol(float(hk_table.cal_val[i]), hk_time, time)
+	return, value
+end
\ No newline at end of file

metis_l1_prep.pro

+pro metis_l1_prep
+
+	; start the log 
+
+	journal,'output/metis_l1_prep_log.txt'
+
+	; some definitions
+
+	metis_datatype = list( $
+		'vl-image', $				; 0
+		'uv-image', $				; 1
+		'pcu-accumul', $			; 2
+		'vl-temp-matrix', $			; 3
+		'uv-temp-matrix', $			; 4
+		'vl-c-ray-mat', $			; 5
+		'uv-c-ray-mat', $			; 6
+		'pcu-event-list', $			; 7
+		'pcu-test-event-list', $	; 8
+		'light-curve' $				; 9
+	)
+
+	metis_obj_size = list( $
+		2048L, $					; 0
+		1024L, $					; 1
+		1024L, $					; 2
+		2048L, $					; 3
+		1024L, $					; 4
+		2048L, $					; 5
+		1024L $						; 6
+	)
+
+	; read the auxiliary file - here we have all the inputs we need
+
+	input = json_parse('input/contents.json', /toarray, /tostruct)
+
+	; prepare the spice kernels
+
+	kernels = [ $
+		input.tls_file_name, $
+		input.tsc_file_name $
+	]
+
+	; load the spice kernels
+
+	cspice_furnsh, kernels
+
+	; read l0 fits structure
+
+	fits_info, input.l0_file_name, n_ext = n_ext, extname = extname, /silent
+
+	; read the primary hdu
+	
+	ext_no = 0
+	image = mrdfits(input.l0_file_name, ext_no, primary_header, /silent)
+
+	; if the data is not an image, read the data binary-table extension
+
+	if fxpar(primary_header, 'NAXIS') eq 0 then begin
+		ext_no += 1
+		data_bin_table = mrdfits(input.l0_file_name, ext_no, data_extension_header, /silent)
+	endif else data_bin_table = !null
+
+	; read the metadata extension
+
+	ext_no += 1
+	metadata_bin_table = mrdfits(input.l0_file_name, ext_no, metadata_extension_header, /silent)
+
+	; identify the data product and its nominal size
+
+	datatype = fxpar(metadata_extension_header, 'DATATYPE')
+
+	if datatype le 6 then begin
+		naxis = 2
+		naxis1 = metis_obj_size[datatype]
+		naxis2 = naxis1
+		bitpix = 16
+		data_size = naxis1 * naxis2
+		data_volume = (data_size * bitpix) / 8.
+	endif else begin
+		naxis = 0
+		bitpix = 8
+		data_size = 0
+		data_volume = 0
+	endelse
+
+	; if the data product is an accumulation matrix, look for the accumulation vector extension and read it if it exists
+
+	if datatype eq 2 then begin
+		if fxpar(metadata_extension_header, 'M') eq 256 then begin
+			ext_no += 1
+			accumul_vector = mrdfits(input.l0_file_name, ext_no, vector_extension_header, /silent)
+		endif else accumul_vector = !null
+	endif
+
+	; NOTE - the radialization extension is ignored
+
+	; read the planning data
+
+	planning_data = json_parse(input.planning_file_name, /toarray, /tostruct)
+
+	; definitions for the primary header
+	; version of the fits file
+
+	version = string(input.l1_version + 1, format = '(I02)')
+
+	; creation and acquisition times in utc
+
+	date = date_conv(systime(/julian, /utc), 'FITS') + 'Z'
+
+	obt_beg = fxpar(primary_header, 'OBT_BEG')
+	obt_end = fxpar(primary_header, 'OBT_END')
+	obt_avg = (obt_beg + obt_end) / 2.0D
+
+	date_beg = solo_obt2utc(decode_obt(obt_beg, /from_double))
+	date_end = solo_obt2utc(decode_obt(obt_end, /from_double))
+	date_avg = solo_obt2utc(decode_obt(obt_avg, /from_double))
+
+	date_beg_string = strreplace(strreplace(strmid(date_beg, 0, 19), '-', ''), ':', '')
+
+	; name of the fits file
+	
+	file_name_fields = strsplit(fxpar(primary_header, 'FILENAME'), '_', /extract)
+	file_name = 'solo_l1_' + file_name_fields[2] + '_' + date_beg_string + '_v' + version + '.fits'
+	out_file_name = 'output/' + file_name
+
+	; instrument keywords
+	; TODO - complete with filter information
+
+	if datatype eq 0 or datatype eq 3 or datatype eq 5 or datatype eq 9 then begin
+		filter = 'VL'
+		wavelnth = 0.0
+		wavemin = 0.0
+		wavemax = 0.0
+		waveband = 'Visible light 580-640 nm'
+	endif else begin
+		filter = 'UV'
+		wavelnth = 1215.67
+		wavemin = 0.0
+		wavemax = 0.0
+		waveband = 'H I Lyman-alpha 121.6 nm'
+	endelse
+	detector = filter + 'DA'
+	telescope = 'SOLO/Metis/' + detector
+
+	; campaign keywords
+
+	obs_id_fields = strsplit(planning_data.obs_id, '_', /extract)
+	soop_type = obs_id_fields[2]
+	obs_type = obs_id_fields[4]
+
+	; build the fits file extensions
+
+	; join the metadata extension header to the primary header removing unwanted keywords
+
+	i = where(strmid(primary_header, 0, 8) eq 'COMP_RAT')
+	j = where(strmid(metadata_extension_header, 0, 8) eq 'EXTNAME ')
+	n = n_elements(metadata_extension_header)
+
+	primary_header = [ $
+		primary_header[0 : i], $
+		metadata_extension_header[j + 1 : n - 2], $
+		primary_header[i + 1 : *] $
+	]
+	
+	; adjust the primary header (it is almost the same for all data product types)
+
+	fxaddpar, primary_header, 'FILENAME', file_name, 'FITS filename'
+	fxaddpar, primary_header, 'PARENT', file_basename(input.l0_file_name), 'Name of the parent file that got processed to the current one', before = 'APID'
+	fxaddpar, primary_header, 'DATE', date, 'Date and time of FITS file creation', before = 'OBT_BEG'
+	fxaddpar, primary_header, 'DATE-OBS', date_beg, 'Same as DATE-BEG', before = 'OBT_BEG'
+	fxaddpar, primary_header, 'DATE-BEG', date_beg, 'Start time of observation', before = 'OBT_BEG'
+	fxaddpar, primary_header, 'DATE-AVG', date_avg, 'Average time of observation', before = 'OBT_BEG'
+	fxaddpar, primary_header, 'DATE-END', date_end, 'End time of observation', before = 'OBT_BEG'
+	fxaddpar, primary_header, 'TIMESYS', 'UTC', 'System used for time keywords', before = 'OBT_BEG'
+	fxaddpar, primary_header, 'TIMRDER', 0.0, 'Estimated random error in time values', before = 'OBT_BEG'
+	fxaddpar, primary_header, 'TIMSYER', 0.0, 'Estimated systematic error in time values', before = 'OBT_BEG'
+	fxaddpar, primary_header, 'LEVEL', 'L1', 'Data processing level'
+	fxaddpar, primary_header, 'CREATOR', 'metis_l1_prep.pro', 'FITS creation software'
+	fxaddpar, primary_header, 'OBSRVTRY', 'Solar Orbiter', 'Satellite name', before = 'INSTRUME'
+	fxaddpar, primary_header, 'TELESCOP', telescope, 'Telescope that took the measurement', before = 'INSTRUME'
+	fxaddpar, primary_header, 'DETECTOR', detector, 'Subunit/sensor', before = 'DATAMIN'
+	fxaddpar, primary_header, 'OBJECT', 'TBD', 'The use of the keyword OBJECT is [TBD]', before = 'DATAMIN'
+	fxaddpar, primary_header, 'OBS_MODE', planning_data.obs_mode, 'Observation mode or study that has been used to acquire this image', before = 'DATAMIN'
+	fxaddpar, primary_header, 'OBS_TYPE', obs_type, 'Subfield of OBS_ID that only contains an encoded version of OBS_MODE', before = 'DATAMIN'
+	fxaddpar, primary_header, 'FILTER', filter, 'Filter used to acquire this image', before = 'DATAMIN'
+	fxaddpar, primary_header, 'WAVELNTH', wavelnth, 'Characteristic wavelength at which the observation was taken', before = 'DATAMIN'
+	fxaddpar, primary_header, 'WAVEMIN', wavemin, 'The shortest wavelength at which the net (approximate) response function becomes 0.05 times the maximum response. ', before = 'DATAMIN'
+	fxaddpar, primary_header, 'WAVEMAX', wavemax, 'The longest wavelength at which the net (approximate) response function becomes 0.05 times the maximum response', before = 'DATAMIN'
+	fxaddpar, primary_header, 'WAVEBAND', waveband, 'Description of the wavelength band', before = 'DATAMIN'
+	fxaddpar, primary_header, 'XPOSURE', 0, 'Total effective exposure time of the observation, in seconds', before = 'DATAMIN'
+	fxaddpar, primary_header, 'NSUMEXP', 0, 'Number of images summed together to form the observation', before = 'DATAMIN'
+	fxaddpar, primary_header, 'TELAPSE', 0, 'Total elapsed time between the beginning and end of the complete observation in seconds, including any dead times between exposures', before = 'DATAMIN'
+	fxaddpar, primary_header, 'SOOPNAME', planning_data.soop_name, 'SOOP(s) that this observation belongs to', before = 'DATAMIN'
+	fxaddpar, primary_header, 'SOOPTYPE', soop_type, before = 'DATAMIN'
+	fxaddpar, primary_header, 'OBS_ID', planning_data.obs_id, 'Unique identifier for the observation that is associated with the data acquisition', before = 'DATAMIN' 
+	fxaddpar, primary_header, 'TARGET', 'TBD', 'Taget as defined in the SOOP description', before = 'DATAMIN'
+	fxaddpar, primary_header, 'BSCALE', 1.0, before = 'DATAMIN'
+	fxaddpar, primary_header, 'BZERO', 0, before = 'DATAMIN'
+	fxaddpar, primary_header, 'BTYPE', strupcase(metis_datatype[datatype]), 'Science data object type', before = 'DATAMIN'
+	fxaddpar, primary_header, 'BUNIT', 'DN', 'Units of physical value, after application of BSCALE and BZERO', before = 'DATAMIN'
+	if datatype le 6 then begin
+		fxaddpar, primary_header, 'PXBEG1', 1, 'First pixel that has been read out in dimension 1', before = 'COMPRESS'
+		fxaddpar, primary_header, 'PXBEG2', 1, 'First pixel that has been read out in dimension 2', before = 'COMPRESS'
+		fxaddpar, primary_header, 'PXEND1', naxis1, 'Last pixel that has been read out in dimension 1', before = 'COMPRESS'
+		fxaddpar, primary_header, 'PXEND2', naxis2, 'Last pixel that has been read out in dimension 2', before = 'COMPRESS'
+	endif
+	fxaddpar, primary_header, 'NBIN1', 1, 'Binning factor in the dimension 1', before = 'COMPRESS'
+	fxaddpar, primary_header, 'NBIN2', 1, 'Binning factor in the dimension 2', before = 'COMPRESS'
+	fxaddpar, primary_header, 'NBIN', 1, 'Product of all NBIN values above', before = 'COMPRESS'
+	fxaddpar, primary_header, 'IDB_VERS', input.idb_version, '', before = 'HDR_VERS'
+	fxaddpar, primary_header, 'INFO_URL', 'http://metis.oato.inaf.it', 'Link to more information on the instrument data', before = 'HISTORY'
+
+	; add checksum and datasum to the fits header
+	; WARN - should this be done at the end? i don't know
+
+	fits_add_checksum, primary_header, image
+
+	; add keywords for file history
+
+	fxaddpar, primary_header, 'HISTORY', 'Third history line'
+
+	; remove unused keywords
+
+	sxdelpar, primary_header, 'CAD_BEG'
+	sxdelpar, primary_header, 'CAD_END'
+	sxdelpar, primary_header, 'WIDTH'
+	sxdelpar, primary_header, 'HEIGHT'
+	sxdelpar, primary_header, 'NB_IMG'
+	sxdelpar, primary_header, 'N'
+	sxdelpar, primary_header, 'X_SIZE'
+	sxdelpar, primary_header, 'Y_SIZE'
+	sxdelpar, primary_header, 'Z_SIZE'
+	sxdelpar, primary_header, 'P_BANDS'
+	sxdelpar, primary_header, 'N_BANDS'
+	sxdelpar, primary_header, 'ORIG_X'
+	sxdelpar, primary_header, 'ORIG_Y'
+	sxdelpar, primary_header, 'FIRSTROW'
+
+	mwrfits, image, out_file_name, primary_header, /no_comment, /create, /silent
+
+	; if applicable, save the data binary-table extension as it is
+
+	if isa(data_bin_table) then mwrfits, data_bin_table, 'output/' + file_name, data_extension_header, /no_comment, /silent
+
+	; read the house-keeping telemetry
+
+	hk_bin_table = json_parse(input.hk_file_name, /toarray, /tostruct)
+
+	; replace raw values with calibrated values in the primary header
+	; WARN - must be done here?
+
+	if datatype eq 0 or datatype eq 3 or datatype eq 5 then begin
+		fxaddpar, primary_header, 'DAC1POL1', interpol_param(hk_bin_table, 'NIT0E061', obt_avg)
+		fxaddpar, primary_header, 'DAC2POL1', interpol_param(hk_bin_table, 'NIT0E062', obt_avg)
+		fxaddpar, primary_header, 'DAC1POL2', interpol_param(hk_bin_table, 'NIT0E064', obt_avg)
+		fxaddpar, primary_header, 'DAC2POL2', interpol_param(hk_bin_table, 'NIT0E065', obt_avg)
+		fxaddpar, primary_header, 'DAC1POL3', interpol_param(hk_bin_table, 'NIT0E067', obt_avg)
+		fxaddpar, primary_header, 'DAC2POL3', interpol_param(hk_bin_table, 'NIT0E068', obt_avg)
+		fxaddpar, primary_header, 'DAC1POL4', interpol_param(hk_bin_table, 'NIT0E06A', obt_avg)
+		fxaddpar, primary_header, 'DAC2POL4', interpol_param(hk_bin_table, 'NIT0E06B', obt_avg)
+		fxaddpar, primary_header, 'TSENSOR', interpol_param(hk_bin_table, 'NIT0E0E0', obt_avg)
+		fxaddpar, primary_header, 'PMPTEMP', interpol_param(hk_bin_table, 'NIT0L00D', obt_avg)
+	endif
+
+	
+	if datatype eq 1 or datatype eq 4 or datatype eq 6 then begin
+		fxaddpar, primary_header, 'HVU_SCR', interpol_param(hk_bin_table, 'NIT0E070', obt_avg)
+		fxaddpar, primary_header, 'HVU_MCP', interpol_param(hk_bin_table, 'NIT0E071', obt_avg)
+		fxaddpar, primary_header, 'TSENSOR', interpol_param(hk_bin_table, 'NIT0E050', obt_avg)
+		
+	endif
+
+	; build the telemetry extension
+
+	hk_extension_header = !null
+	fxaddpar, hk_extension_header, 'XTENSION', 'BINTABLE', 'Extension type'
+	fxaddpar, hk_extension_header, 'BITPIX', 8, 'Number of bits per data pixel'
+	fxaddpar, hk_extension_header, 'NAXIS', 2, 'Number of data axes'
+	fxaddpar, hk_extension_header, 'NAXIS1', 0, 'Number of pixels along axis 1'
+	fxaddpar, hk_extension_header, 'NAXIS2', 0, 'Number of pixels along axis 2'
+	fxaddpar, hk_extension_header, 'PCOUNT', 0, 'Parameter count'
+	fxaddpar, hk_extension_header, 'GCOUNT', 1, 'Group count'
+	fxaddpar, hk_extension_header, 'EXTNAME', 'House-keeping', 'Extension name'
+
+	mwrfits, hk_bin_table, out_file_name, hk_extension_header, /no_comment, /silent
+
+	; write the auxiliary information file
+
+	output = { $
+		file_name: out_file_name, $
+		l0_file_name: input.l0_file_name, $
+		log_file_name: '' $
+	}
+	
+	openw, unit, 'output/contents.json'
+	printf, unit, output, /implied_print
+	free_lun, unit
+	
+	; json_write, output, 'output/contents.json'
+
+	; unload the spice kernels
+
+	cspice_unload, kernels
+
+	; close the log
+
+	journal
+end
\ No newline at end of file

solo_obt2utc.pro

+function solo_obt2utc, obt
+	solo_id = -144
+	cspice_scs2e, solo_id, obt, et
+	cspice_et2utc, et, 'ISOC', 2, utc
+	return, utc + 'Z'
+end
\ No newline at end of file

strreplace.pro

+function strreplace, input, substring, replacement
+	return, strjoin(strsplit(input, substring, /extract, /preserve), replacement)
+end