Created
October 10, 2022 21:53
-
-
Save cosmonaut/af9839c7b6492eceaf3ac9cc96238962 to your computer and use it in GitHub Desktop.
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| # Author: Nicholas Nell | |
| # email: [email protected] | |
| # | |
| # MetroPro DAT file extraction script -- this program extracts the | |
| # basic header and intensity/phase information contained in a MetroPro | |
| # DAT binary file. | |
| # | |
| # File information from "MetroPro Reference Guide" section 12 for the | |
| # "MetroPro Binary Data File Format" | |
| # | |
| # This script outputs images of the intensity data and phase data if | |
| # they exist (with suffixes .intsty.png and .phase.png attached to the | |
| # DAT filename). It also saves .npy files of the intensity and phase | |
| # data with suffixes .intsty.npy and .phase.npy appended to the | |
| # original DATA filename. | |
| import sys | |
| import numpy as np | |
| import struct | |
| import matplotlib.pyplot as mp | |
| def main(fname): | |
| with open(fname, mode = 'rb') as f: | |
| data = f.read() | |
| # "Nominal" header size in bytes | |
| hdr_size = 834 | |
| if (len(data) < 834): | |
| print("File too small to be MetroPro DAT format") | |
| exit() | |
| # Read first 4 bytes to determine "magic number" | |
| magic_number = data[0:4] | |
| print("MAGIC NUMBER: 0x{0:s}\n".format(magic_number.hex())) | |
| if (magic_number == bytearray([0x88, 0x1B, 0x03, 0x71])): | |
| print("Detected header format 3") | |
| hdr_size = 4096 | |
| elif (magic_number == bytearray([0x88, 0x1B, 0x03, 0x70])): | |
| print("Detected header format 2") | |
| elif (magic_number == bytearray([0x88, 0x1B, 0x03, 0x6f])): | |
| print("Detected header format 1") | |
| else: | |
| print("No valid magic number detected, may not be dat file format") | |
| exit() | |
| hdr_format = struct.unpack('>h', data[4:6])[0] | |
| print("HEADER FORMAT: {0:d}".format(hdr_format)) | |
| hsize = struct.unpack('>l', data[6:10])[0] | |
| print("HEADER SIZE: {0:d}".format(hsize)) | |
| # Intensity data matrix header variables | |
| ac_org_x = struct.unpack('>h', data[48:50])[0] | |
| ac_org_y = struct.unpack('>h', data[50:52])[0] | |
| ac_width = struct.unpack('>h', data[52:54])[0] | |
| ac_height = struct.unpack('>h', data[54:56])[0] | |
| ac_n_buckets = struct.unpack('>h', data[56:58])[0] | |
| ac_range = struct.unpack('>h', data[58:60])[0] | |
| ac_n_bytes = struct.unpack('>l', data[60:64])[0] | |
| #print(ac_n_bytes) | |
| print("") | |
| print("INTENSITY DATA MATRIX HEADER DATA") | |
| print("X/Y ORIGIN: {0:d}, {1:d}".format(ac_org_x, ac_org_y)) | |
| print("WIDTH: {0:d}".format(ac_width)) | |
| print("HEIGHT: {0:d}".format(ac_height)) | |
| print("N BUCKETS: {0:d}".format(ac_n_buckets)) | |
| print("RANGE: {0:d}".format(ac_range)) | |
| print("AC TOTAL BYTES: {0:d}".format(ac_n_bytes)) | |
| # If the intensity data matrix is present then it begins | |
| # immediately after the header | |
| if (ac_n_bytes > 0): | |
| imtrx = np.frombuffer(data[hdr_size:hdr_size + ac_n_bytes], dtype = '>u2').copy() | |
| imtrx = imtrx.reshape(ac_height, ac_width) | |
| # Supress 0xffff (masked) bytes | |
| ind = imtrx > ac_range | |
| imtrx[ind] = 0 | |
| # Generate image | |
| mp.imshow(imtrx) | |
| mp.colorbar() | |
| mp.title("Intensity Map") | |
| mp.savefig(fname + ".intsty.png") | |
| # Save numpy file | |
| np.save(fname + ".intsty", imtrx) | |
| mp.clf() | |
| # Phase data matrix header variables | |
| cn_org_x = struct.unpack('>h', data[64:66])[0] | |
| cn_org_y = struct.unpack('>h', data[66:68])[0] | |
| cn_width = struct.unpack('>h', data[68:70])[0] | |
| cn_height = struct.unpack('>h', data[70:72])[0] | |
| cn_n_bytes = struct.unpack('>l', data[72:76])[0] | |
| intf_scale_factor = struct.unpack('>f', data[164:168])[0] | |
| wavelength_in = struct.unpack('>f', data[168:172])[0] | |
| obliquity_factor = struct.unpack('>f', data[176:180])[0] | |
| phase_res = struct.unpack('>h', data[218:220])[0] | |
| print("") | |
| print("PHASE MATRIX HEADER DATA") | |
| print("X/Y ORIGIN: {0:d}, {1:d}".format(cn_org_x, cn_org_y)) | |
| print("WIDTH {0:d}".format(cn_width)) | |
| print("HEIGHT {0:d}".format(cn_height)) | |
| print("INTF_SCALE_FACTOR: {0:f}".format(intf_scale_factor)) | |
| print("WAVELENGTH IN: {0:.6g} (m)".format(wavelength_in)) | |
| print("OBLIQUITY FACTOR: {0:f}".format(obliquity_factor)) | |
| print("PHASE RES: {0:d}".format(phase_res)) | |
| if (hdr_format == 1): | |
| if (phase_res == 0): | |
| R = 4096 | |
| elif (phase_res == 1): | |
| R = 32768 | |
| else: | |
| print("INVALID PHASE RES") | |
| exit() | |
| elif ((hdr_format == 2) or (hdr_format == 3)): | |
| if (phase_res == 0): | |
| R = 4096 | |
| elif (phase_res == 1): | |
| R = 32768 | |
| elif (phase_res == 2): | |
| R = 131072 | |
| else: | |
| print("INVALID PHASE RES") | |
| exit() | |
| else: | |
| print("Invalid header format.") | |
| exit() | |
| print("PHASE RES ELEMENTS (R): {0:d}".format(R)) | |
| S = intf_scale_factor | |
| O = obliquity_factor | |
| # If phase data matrix is present process it. | |
| if (cn_n_bytes > 0): | |
| pmtrx = np.frombuffer(data[hdr_size + ac_n_bytes:hdr_size + ac_n_bytes + cn_n_bytes], dtype = '>i4').copy() | |
| pmtrx = pmtrx.reshape(cn_height, cn_width) | |
| ind = pmtrx >= 0x7ffffff8 | |
| pmtrx[ind] = 0 | |
| # This is now the pmtrx in units of "zygos" | |
| # Convert to units of waves | |
| pmtrx = pmtrx*((S*O)/R) | |
| # It is also possible to use the formula (W*S*O)/R to have the | |
| # units in meters where W is wavelength is given in the | |
| # "wavelength_in" header variable in units of meters. | |
| # Generate image | |
| mp.imshow(pmtrx) | |
| mp.colorbar() | |
| mp.title("Phase Map (waves)") | |
| mp.savefig(fname + ".phase.png") | |
| # Save numpy file | |
| np.save(fname + ".phase", imtrx) | |
| # User can continue to use imtrx and pmtrx here as desired. | |
| if __name__ == '__main__': | |
| if (len(sys.argv) < 1): | |
| print("Usage: metropro_extract.py <filename>") | |
| exit() | |
| main(sys.argv[1]) |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment