Commit b9b676d1 authored by Freddie Akeroyd's avatar Freddie Akeroyd
Browse files

Add first files ... need major tidy ups so not to be taken too seriously

as yet
Refs #52.
parent f7aba55d
/*
C
C Compression of 32bit integer data into byte-relative format
C
C Each integer is stored relative to the previous value in byte form.The first
C is relative to zero. This allows for numbers to be within + or - 127 of
C the previous value. Where a 32bit integer cannot be expressed in this way a
C special byte code is used (-128) and the full 32bit integer stored elsewhere.
C The final space used is (NIN-1)/4 +1 + NEXTRA longwords, where NEXTRA is the
C number of extra longwords used in giving absolute values.
C
*/
#include <stdexcept>
#include <iostream>
#include "byte_rel_comp.h"
#define LARGE_NUMBER 1073741824
#define FAILURE 1
#define SUCCESS 0
int byte_rel_comp(int* data_in, int n_in, char* data_out, int max_out, int& n_out)
{
int i, icurrent, irel;
union
{
int i;
char c[4];
} byte_pack;
if (n_in <= 0)
{
throw std::runtime_error("byte rel comp error: nin <= 0");
}
if (max_out <= n_in)
{
throw std::runtime_error("byte rel comp error: nin <= 0");
}
n_out = 0;
icurrent = 0;
for(i=0; i<n_in; i++)
{
// Trap out ridiculously large numbers. They could cause problems in subtraction,
// so force them to be stored absolutely
if ( (data_in[i] > LARGE_NUMBER) || (data_in[i] < -LARGE_NUMBER) || (icurrent > LARGE_NUMBER) || (icurrent < -LARGE_NUMBER) )
{
irel = 128; // Force absolute mode
}
else
{
irel = data_in[i] - icurrent; // Calc relative offset
}
// If small put it in a byte
if ( (irel <= 127) && (irel >= -127) )
{
if (n_out > max_out){
throw std::runtime_error("byte rel comp error: nin <= 0");
}
data_out[n_out] = irel; // Pack relative byte
n_out++;
}
else
{
// Otherwise put marker in byte followed by packed 32bit integer
if (n_out + 4 >= max_out)
{
throw std::runtime_error("byte rel comp error: nin <= 0");
}
data_out[n_out] = -128; // pack marker
byte_pack.i = data_in[i];
data_out[n_out+1] = byte_pack.c[0];
data_out[n_out+2] = byte_pack.c[1];
data_out[n_out+3] = byte_pack.c[2];
data_out[n_out+4] = byte_pack.c[3];
n_out += 5;
}
icurrent = data_in[i];
}
return SUCCESS;
}
/*
C
C Expansion of byte-relative format into 32bit format
C
C Each integer is stored relative to the previous value in byte form.The first
C is relative to zero. This allows for numbers to be within + or - 127 of
C the previous value. Where a 32bit integer cannot be expressed in this way a
C special byte code is used (-128) and the full 32bit integer stored elsewhere.
C The final space used is (NIN-1)/4 +1 + NEXTRA longwords, where NEXTRA is the
C number of extra longwords used in giving absolute values.
C
C
C Type definitions
C Passed parameters
INTEGER NIN
BYTE INDATA(NIN)
C Array of NIN compressed bytes
INTEGER NFROM
C pass back data from this 32bit word onwards
C
INTEGER NOUT
C Number of 32bit words expected
INTEGER OUTDATA(NOUT)
C outgoing expanded data
INTEGER ISTATUS
C Status return
C =1 no problems!
C =3 NOUT .lt.NIN/5
C =2 NIN .le.0
C =4 NOUT .gt.NIN
C =6 number of channels lt NOUT
*/
//n_from is zero based
int byte_rel_expn(char* data_in, int n_in, int n_from, int* data_out, int n_out)
{
int i, j;
union
{
int i;
char c[4];
} byte_pack;
// First check no slip-ups in the input parameters
if (n_in <= 0)
{
throw std::runtime_error("byte rel comp error: nin <= 0");
}
if (n_out + n_from > n_in)
{
throw std::runtime_error("byte rel comp error: nin <= 0");
}
// Set initial absolute value to zero and channel counter to zero
byte_pack.i = 0;
j = 0;
// Loop over all expected 32bit integers
for(i=0; i< n_from + n_out; i++)
{
if (j >= n_in)
{
throw std::runtime_error("byte rel comp error: nin <= 0");
}
// if number is contained in a byte
if (data_in[j] != -128)
{
// add in offset to base
byte_pack.i += data_in[j];
j++;
}
else
{
// Else skip marker and pick up new absolute value
// check there are enough bytes
if (j+4 >= n_in)
{
throw std::runtime_error("byte rel comp error: nin <= 0");
}
// unpack 4 bytes
byte_pack.c[0] = data_in[j+1];
byte_pack.c[1] = data_in[j+2];
byte_pack.c[2] = data_in[j+3];
byte_pack.c[3] = data_in[j+4];
// CALL VAX_TO_LOCAL_INTS(ITEMP, 1)
j += 5;
}
// update current value
if (i >= n_from)
{
data_out[i-n_from] = byte_pack.i;
}
}
// expansion OK, but excessive number of bytes given to the routine
if (n_out < n_in / 5)
{
std::cerr << "byte rel expn: excessive bytes" << std::endl;
}
return SUCCESS;
}
int byte_rel_comp(int* data_in, int n_in, char* data_out, int max_out, int& n_out);
int byte_rel_expn(char* data_in, int n_in, int n_from, int* data_out, int n_out);
This diff is collapsed.
#ifndef ISISRAW_H
#define ISISRAW_H
struct ISISCRPT_STRUCT;
#include "item_struct.h"
struct HDR_STRUCT
{
char inst_abrv[3]; // instrument abbreviated name
char hd_run[5]; // run number
char hd_user[20]; // user name
char hd_title[24]; // short title
char hd_date[12]; // start date
char hd_time[8]; // start time
char hd_dur[8]; // run duration (uA.hour)
HDR_STRUCT() { memset(this, ' ', sizeof(HDR_STRUCT)); }
}; // 80 bytes
// address offsets to various sections
struct ADD_STRUCT
{
int ad_run; // 32 (1+offset_of(ver2))
int ad_inst; // ad_run + 94
int ad_se; // ad_inst + 70 + 2*nmon+(5+nuse)*ndet
int ad_dae; // ad_se + 66 + nse*32
int ad_tcb; // ad_dae + 65 + 5*ndet
int ad_user; // ad_tcb + 288 + (ntc1 + 1)
int ad_data; // ad_user + 2 + ulen
int ad_log; // ad_data + 33 +
int ad_end; // 1+end of file
ADD_STRUCT() { memset(this, 0, sizeof(ADD_STRUCT)); }
}; // 9*4 bytes
struct USER_STRUCT
{
char r_user[20]; // name
char r_daytel[20];
char r_daytel2[20];
char r_night[20];
char r_instit[20]; // institute
char unused[3][20]; // to pad out to 8*20 bytes
USER_STRUCT() { memset(this, ' ', sizeof(USER_STRUCT)); }
}; // 8*20 bytes
struct RPB_STRUCT
{
int r_dur; // actual run duration
int r_durunits; // scaler for above (1=seconds)
int r_dur_freq; // test interval for above (seconds)
int r_dmp; // dump interval
int r_dmp_units; // scaler for above
int r_dmp_freq; // test interval for above
int r_freq; // 2**k where source frequency = 50 / 2**k
float r_gd_prtn_chrg; // good proton charge (uA.hour)
float r_tot_prtn_chrg; // total proton charge (uA.hour)
int r_goodfrm; // good frames
int r_rawfrm; // raw frames
int r_dur_wanted; // requested run duration (units as for "duration" above)
int r_dur_secs; // actual run duration in seconds
int r_mon_sum1; // monitor sum 1
int r_mon_sum2; // monitor sum 2
int r_mon_sum3; // monitor sum 3
char r_enddate[12]; // format DD-MMM-YYYY
char r_endtime[8]; // format HH-MM-SS
int r_prop; // RB (proposal) number
int spare[10]; // to pad out to 32*4 bytes
RPB_STRUCT() { memset(this, 0, sizeof(RPB_STRUCT)); }
}; // 32*4 bytes
struct IVPB_STRUCT
{
float i_chfreq; // frequency chopper 1 (Hz)
float freq_c2; // frequency chopper 2 (Hz)
float freq_c3; // frequency chopper 3 (Hz)
int delay_c1; // delay chopper 1 (us)
int delay_c2; // delay chopper 2 (us)
int delay_c3; // delay chopper 3 (us)
int delay_error_c1; // max delay error chopper 1 (us)
int delay_error_c2; // max delay error chopper 2 (us)
int delay_error_c3; // max delay error chopper 3 (us)
float i_chopsiz;
int aperture_c2;
int aperture_c3;
int status_c1; // status c1 (run,stopped,stop open)
int status_c2; // status c2 (run,stopped,stop open)
int status_c3; // status c3 (run,stopped,stop open)
int i_mainshut; // main shutter open = 1
int i_thermshut; // thermal shutter open = 1
float i_xsect; // beam aperture horizontal (mm)
float i_ysect; // beam aperture vertical (mm)
int i_posn; // scattering position (1 or 2, for HRPD)
int i_mod; // moderator type code
int i_vacuum; // detector_tank_vacuum 1 = vacuum on
float i_l1; // L1 scattering length
int i_rfreq; // rotor_frequency (HET)
float i_renergy; // (HET)
float i_rphase; // (HET)
int i_rslit; // slit package (0="unknown",1="L",2="Med",3="Hi") HET
int i_slowchop; // (1=on, 0=off) HET
float i_xcen; // LOQ x centre
float i_ycen; // LOQ y centre
int i_bestop; // beam_stop LOQ
float i_radbest; // beam_stop_radius LOQ
float i_sddist; // source to detector distance (LOQ)
float i_foeang; // foe angle LOQ
float i_aofi; // angle of incidence (CRISP)
int spare[29]; // to pad out to 64*4 bytes
IVPB_STRUCT() { memset(this, 0, sizeof(IVPB_STRUCT)); }
}; // 64*4 bytes
struct SPB_STRUCT
{
int e_posn; // sample changer position
int e_type; // sample type (1=sample+can,2=empty can)
int e_geom; // sample geometry
float e_thick; // sample thickness normal to sample (mm)
float e_height; // sample height (mm)
float e_width; // sample width (mm)
float e_omega; // omega sample angle (degrees)
float e_chi; // chi sample angle (degrees)
float e_phi; // phi sample angle (degrees)
float e_scatt; // scattering geometry (1=trans, 2 =reflect)
float e_xscatt; // sample coherent scattering cross section (barn)
float samp_cs_inc;
float samp_cs_abs;
float e_dens; // sample number density (atoms.A-3)
float e_canthick; // can wall thickness (mm)
float e_canxsect; // can coherent scattering cross section (barn)
float can_cs_inc;
float can_cs_abs;
float can_nd; // can number density (atoms.A-3)
char e_name[40]; // sample name of chemical formula
int e_equip;
int e_eqname;
int spare[33]; // to bring up to 64*4 bytes
SPB_STRUCT() { memset(this, 0, sizeof(SPB_STRUCT)); }
}; // 64*4 bytes
struct SE_STRUCT
{
char sep_name[8]; // SE block name
int sep_value;
int sep_exponent;
char sep_units[8]; // units of value
int sep_low_trip;
int sep_high_trip;
int sep_cur_val;
int sep_status; // are we in bounds
int sep_control; // controlled parameter (true/false)
int sep_run; // run control parameter (true/false)
int sep_log; // logged parameter (true/false)
float sep_stable; // units per sec
float sep_period; // monitor repeat period
int sep_cam_addr; // camac location N
int sep_cam_sub; // camac location A
int sep_offset; // CAMAC offset (added to value)
int sep_cam_rgrp; // camac register group (1 or 2)
int sep_pre_proc; // pre process routine number
int sep_cam_vals[12]; // camac values
SE_STRUCT() { memset(this, 0, sizeof(SE_STRUCT)); }
}; // 32*4 bytes
struct DAEP_STRUCT
{
int a_pars; // Word length in bulk store memory
int mem_size; // Length of bulk store memory (bytes)**A
int a_minppp; // PPP minimum value ***B
int ppp_good_high; // good PPP total (high 32 bits) ***B
int ppp_good_low; // good PPP total (low 32 bits) ***B
int ppp_raw_high; // raw PPP total (high 32 bits) ***B
int ppp_raw_low; // raw PPP total (low 32 bits) ***B
int neut_good_high; // good ext. neut tot (high 32bits)***B
int neut_good_low; // good ext. neut tot (low 32 bits)***B
int neut_raw_high; // raw ext. neut tot (high 32 bits)***B
int neut_raw_low; // raw ext. neut tot (low 32 bits)***B
int neut_gate_t1; // ext. neutron gate (t1) (s) ***B
int neut_gate_t2; // ext. neutron gate (t2) (s) ***B
int mon1_detector; // detector for MON 1 (12 bits) ***B
int mon1_module; // module for MON 1 ( 4 bits) ***B
int mon1_crate; // crate for MON 1 ( 4 bits) ***B
int mon1_mask; // mask for MON 1 (c4:m4:d12) ***B
int mon2_detector; // detector for MON 2 (12 bits) ***B
int mon2_module; // module for MON 2 ( 4 bits) ***B
int mon2_crate; // crate for MON 2 ( 4 bits) ***B
int mon2_mask; // mask for MON 2 (c4:m4:d12) ***B
int events_good_high; // total GOOD EVENTS (high 32 bits)***B
int events_good_low; // total GOOD EVENTS (low 32 bits)***B
int a_delay; // frame synch delay (4s steps) ***B
int a_sync; // frm snch origin(0:none/1:ext/2:int)***B
int a_smp; // Secondary Master Pulse (0:en,1:dis)
int ext_vetos[3]; // External vetoes 0,1,2 (0 dis,1 en)
// extra bits for new PC DAE
int n_tr_shift; // set to number of shifted time regimes
int tr_shift[3]; // set to shift value (us) of each TR if (using_tr_shift > 0)
int spare[31]; // to 64*4 bytes
DAEP_STRUCT() { memset(this, 0, sizeof(DAEP_STRUCT)); }
}; // 64*4 bytes
struct DHDR_STRUCT
{
int d_comp; // compression type (0=none, 1 = byte relative)
int reserved; // unused
int d_offset; // offset to spectrum descriptor array
float d_crdata; // compression ratio for data
float d_crfile; // compression ratio for whole file
int d_exp_filesize; // equivalent version 1 filesize
int unused[26]; // to bring size to 32*4 bytes
DHDR_STRUCT() { memset(this, 0, sizeof(DHDR_STRUCT)); d_comp = 1; d_offset = 1+32; }
}; // 32*4 bytes
struct DDES_STRUCT
{
int nwords; // number of compressed words in spectrum
int offset; // offset to compressed spectrum
DDES_STRUCT() : nwords(0), offset(0) {}
}; // 2*4 bytes
struct LOG_LINE
{
int len; // real length of data
char* data; // padded to multiple of 4 bytes
LOG_LINE() : len(0), data(0) {}
};
struct LOG_STRUCT
{
int ver; // = 2
int nlines;
LOG_LINE* lines; // size nlines
LOG_STRUCT() : ver(2), nlines(0), lines(0) {}
};
class ISISRAW
{
private:
ISISCRPT_STRUCT* m_crpt;
// these are used for caching on updates
int m_ntc1;
int m_nsp1;
int m_nper;
item_struct<char> m_char_items;
item_struct<float> m_real_items;
item_struct<int> m_int_items;
int addItems();
public:
// section 1
HDR_STRUCT hdr; // header block (80 bytes)
int frmt_ver_no; // format version number VER1 (=2)
struct ADD_STRUCT add; // 9*4 bytes
int data_format; // data section format (0 = by TC, 1 = by spectrum)
// section 2
int ver2; // run section version number VER2 (=1)
int r_number; // run number
char r_title[80]; // run title
USER_STRUCT user; // user information (8*20 bytes)
RPB_STRUCT rpb; // run parameter block (32*4 bytes)
// section 3
int ver3; // instrument section version number (=2)
char i_inst[8]; // instrument name
IVPB_STRUCT ivpb; // instrument parameter block (64*4 bytes)
int i_det; // number of detectors NDET
int i_mon; // number of monitors NMON
int i_use; // number of user defined UTn tables NUSE
// I_TABLES is address of MDET
int* mdet; // detector number for monitors (size NMON)
int* monp; // prescale value for each monitor (size NMON)
int* spec; // spectrum number table (size NDET)
float* delt; // hold off table (size NDET)
float* len2; // L2 table (size NDET)
int* code; // code for UTn tables (size NDET)
float* tthe; // 2theta scattering angle (size NDET)
float* ut; // nuse UT* user tables (total size NUSE*NDET) ut01=phi
// section 4
int ver4; // SE section version number (=2)
SPB_STRUCT spb; // sample parameter block (64*4 bytes)
int e_nse; // number of controlled SEPs NSEP
SE_STRUCT* e_seblock; // NSEP SE parameter blocks (total size NSEP*32*4 bytes)
// section 5
int ver5; // DAE section version number (=2)
DAEP_STRUCT daep; // DAE parameter block (size 64*4 bytes)
// A_TABLES points at CRAT
int* crat; // crate number for each detector (size NDET)
int* modn; // module number for each detector (size NDET)
int* mpos; // module position for each detector (size NDET)
int* timr; // time regime for each detector (size NDET)
int* udet; // user detector number for each detector (size NDET)
// section 6
int ver6; // TCB section version number (=1)
int t_ntrg; // number of time regimes (=1)
int t_nfpp; // number of frames per period
int t_nper; // number of periods
int t_pmap[256]; // period number for each basic period
int t_nsp1; // number of spectra in time regime 1
int t_ntc1; // number of time channels in time regime 1
int t_tcm1[5]; // time channel mode
float t_tcp1[5][4]; // time channel parameters
int t_pre1; // prescale for 32MHz clock
int* t_tcb1; // time channel boundaries in clock pulses (size NTC1+1)
// section 7
// write NCHAN = NTC1+1 time channel boundaries
int ver7; // user version number (=1)
int u_len;
float* u_dat; // user defined data (ULEN, max size 400 words)
// section 8
int ver8; // data version number (=2)
DHDR_STRUCT dhdr; // size 32*4 bytes
// D_DATA points at ddes
DDES_STRUCT* ddes; // (NSP1+1)*NPER items, totoal size (NSP1+1)*NPER*2*4 bytes
uint32_t* dat1; // compressed data for (NTC1+1)*(NSP1+1)*NPER values
LOG_STRUCT logsect;
public:
static int vmstime(char* timbuf, int len, time_t time_value);
int size_check();
ISISRAW();
ISISRAW(ISISCRPT_STRUCT* crpt);
int updateFromCRPT();
int ioRAW(FILE* file, bool from_file);
int ioRAW(FILE* file, HDR_STRUCT* s, int len, bool from_file);
int ioRAW(FILE* file, ADD_STRUCT* s, int len, bool from_file);
int ioRAW(FILE* file, USER_STRUCT* s, int len, bool from_file);
int ioRAW(FILE* file, RPB_STRUCT* s, int len, bool from_file);
int ioRAW(FILE* file, IVPB_STRUCT* s, int len, bool from_file);
int ioRAW(FILE* file, SPB_STRUCT* s, int len, bool from_file);
int ioRAW(FILE* file, SE_STRUCT* s, int len, bool from_file);
int ioRAW(FILE* file, DAEP_STRUCT* s, int len, bool from_file);
int ioRAW(FILE* file, DHDR_STRUCT* s, int len, bool from_file);
int ioRAW(FILE* file,
DDES_STRUCT* s, int len, bool from_file);
int ioRAW(FILE* file, LOG_STRUCT* s, int len, bool from_file);
int ioRAW(FILE* file, LOG_LINE* s, int len, bool from_file);
int ioRAW(FILE* file, char* val, int len, bool from_file);
int ioRAW(FILE* file, int* val, int len, bool from_file);
int ioRAW(FILE* file, uint32_t* val, int len, bool from_file);
int ioRAW(FILE* file, float* val, int len, bool from_file);
// allocate
int ioRAW(FILE* file, char** val, int len, bool from_file);
int ioRAW(FILE* file, int** val, int len, bool from_file);
int ioRAW(FILE* file, uint32_t** val, int len, bool from_file);
int ioRAW(FILE* file, float** val, int len, bool from_file);
int ioRAW(FILE* file, SE_STRUCT** s, int len, bool from_file);
int ioRAW(FILE* file, DDES_STRUCT** s, int len, bool from_file);
int ioRAW(FILE* file, LOG_LINE** s, int len, bool from_file);
int readFromFile(const char* filename);
int writeToFile(const char* filename);
int printInfo(std::ostream& os);
};
#endif /* ISISRAW_H */
#include "item_struct.h"
#define FAILURE 1
#define SUCCESS 0
template <typename T>
int item_struct<T>::getItem(const std::string& item_name, T& value)
{
int n;
long spec_no;
// handle case of item_spectrum which means we average the array
// item[ndet] for that particular spectrum
n = item_name.find('_');
if (n != std::string::npos)
{
spec_no = atol(item_name.c_str()+n+1);
return getItem(item_name.substr(0,n), &spec_no, 1, &value);
}
else
{
spec_no = 0;
return getItem(item_name, &spec_no, 0, &value);
}
}
// nspec number of 0 means no spec average
template <typename T>
int item_struct<T>::getItem(const std::string& item_name, long* spec_array, int nspec, T* lVal)
{
int i, j, n;
std::string sitem_name;
std::string tmp_item;
const T* pVal = NULL;
const item_t* item;
item = findItem(item_name, false);
if (item != NULL)
{
for(j=0; j < (nspec == 0 ? 1 : nspec); j++)
{
lVal[j] = *(item->value);
}
return SUCCESS;
}
if (nspec == 0)
{
return FAILURE;
}
item = findItem(item_name, true);
if (item != NULL)
{