Skip to content
Snippets Groups Projects
Unverified Commit 5228b652 authored by acpaquette's avatar acpaquette Committed by GitHub
Browse files

Added themis notebook to crop themis images (#4204)

parent 820a2594
No related branches found
No related tags found
No related merge requests found
%% Cell type:code id: tags:
``` python
import os.path
import datetime
import pvl
import numpy as np
import matplotlib.pyplot as plt
```
%% Cell type:code id: tags:
``` python
themis_file = '/Path/to/themis/raw_data.QUB'
image_file = themis_file
```
%% Cell type:code id: tags:
``` python
header = pvl.load(image_file)
header
```
%% Cell type:code id: tags:
``` python
with open(image_file, 'rb') as f:
# -1 offset obtained from ISIS, record offset likely one based
# so subtract one as we need it to be zero based
image_offset = int((header['^SPECTRAL_QUBE'] - 1) * header['RECORD_BYTES'])
f.seek(image_offset)
b_image_data = f.read()
```
%% Cell type:code id: tags:
``` python
n_bands = header['SPECTRAL_QUBE']['CORE_ITEMS'][header['SPECTRAL_QUBE']['AXIS_NAME'].index('BAND')]
# n_bands = 5
n_lines = header['SPECTRAL_QUBE']['CORE_ITEMS'][header['SPECTRAL_QUBE']['AXIS_NAME'].index('LINE')]
n_lines = 10
if 'LINE_SUFFIX_ITEM_BYTES' in header['SPECTRAL_QUBE'].keys():
n_lines += 2
line_length = header['RECORD_BYTES']
# Add line suffix offset to handle reading two extra lines from each band
if 'LINE_SUFFIX_ITEM_BYTES' in header['SPECTRAL_QUBE'].keys():
bytes_per_band = ((header['SPECTRAL_QUBE']['CORE_ITEMS'][header['SPECTRAL_QUBE']['AXIS_NAME'].index('LINE')]) + header['SPECTRAL_QUBE']['LINE_SUFFIX_ITEM_BYTES']) * line_length
else:
bytes_per_band = (header['SPECTRAL_QUBE']['CORE_ITEMS'][1]) * line_length
```
%% Cell type:code id: tags:
``` python
image_data = []
# Custom numpy data type to handle appropriate DN conversion
# Not necessary but was done for debugging purposes
dt = np.dtype(np.int16)
dt = dt.newbyteorder('>')
for i in range(n_bands):
# image_data.append([])
for j in range(n_lines):
start = (j*line_length) + (bytes_per_band * i)
stop = ((j+1)*line_length) + (bytes_per_band * i)
image_sample = np.frombuffer(b_image_data[start:stop], dtype=dt, count=line_length//2)
image_data.append(image_sample)
image_data = np.array(image_data, dtype=dt)
```
%% Cell type:code id: tags:
``` python
plt.figure(0, figsize=(20, 20))
plt.imshow(image_data)
plt.show()
```
%% Cell type:code id: tags:
``` python
class RealIsisCubeLabelEncoder(pvl.encoder.ISISEncoder):
def encode_time(self, value):
if value.microsecond:
second = u'%02d.%06d' % (value.second, value.microsecond)
else:
second = u'%02d' % value.second
time = u'%02d:%02d:%s' % (value.hour, value.minute, second)
return time
```
%% Cell type:code id: tags:
``` python
image_fn, image_ext = os.path.splitext(image_file)
crop = '_cropped'
themis_image_fn = image_fn + crop + image_ext
themis_image_bn = os.path.basename(themis_image_fn)
grammar = pvl.grammar.ISISGrammar()
grammar.comments+=(("#", "\n"), )
encoder = RealIsisCubeLabelEncoder()
if header['SPECTRAL_QUBE']['DESCRIPTION'] == '':
header['SPECTRAL_QUBE']['DESCRIPTION'] = " "
# Overwrite the number of lines in the label
if 'LINE_SUFFIX_ITEM_BYTES' in header['SPECTRAL_QUBE'].keys():
n_lines -= 2
header['SPECTRAL_QUBE']['CORE_ITEMS'][header['SPECTRAL_QUBE']['AXIS_NAME'].index('LINE')] = n_lines
# Calculate the new offset
# Run this twice as the change to the header will change the headers length.
# Then add two to handle an extra new line character later on
header['^SPECTRAL_QUBE'] = pvl.collections.Units(len(pvl.dumps(header, encoder=encoder, grammar=grammar)), 'BYTES')
header['^SPECTRAL_QUBE'] = pvl.collections.Units(len(pvl.dumps(header, encoder=encoder, grammar=grammar)) + 2, 'BYTES')
```
%% Cell type:code id: tags:
``` python
label_fn, label_ext = os.path.splitext(themis_file)
out_label = label_fn + crop + label_ext
pvl.dump(header, out_label, encoder=encoder, grammar=grammar)
```
%% Cell type:code id: tags:
``` python
with open(themis_image_fn, 'ab+') as f:
b_reduced_image_data = image_data.tobytes()
f.seek(0, 2)
f.write(b'\n')
f.write(b_reduced_image_data)
```
%% Cell type:code id: tags:
``` python
new_header = pvl.load(out_label)
new_header
```
%% Cell type:code id: tags:
``` python
```
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment