// --------------------------------------------------------------------------- 
//
// Copyright, Attribute-Assisted Seismic Processing and Interpretation
//            University of Oklahoma, Norman, OK, USA
//
//            Royalty free use for AASPI sponsors and coinvestigators for use 
//                  in research, exploration, with partners, host governments,
//                  and for provision of processing/interpretation service to 
//                  sponsor clients
//
//            Redistribution, sale, or inclusion of this software in software 
//            products outside the sponsor worksite requires a separate         
//            commmercialization agreement with the University of Oklahoma.
//                              
// Author: J. Timothy Kwiatkowski (jtk@ou.edu)
//         Attribute Assisted Seismic Processing and Interpretation
//         University of Oklahoma
//
// Send errors or bug fixes to kmarfurt@ou.edu
//
// --------------------------------------------------------------------------- 

// --------------------------------------------------------------------------- 
// AASPI (SEP/RSF) to SEGY format converter.
// --------------------------------------------------------------------------- 

#include <iostream>
#include <fstream>
#include <sstream>
#include <cstdlib>
#include <set>
#include <cmath>
#include <algorithm>
#include "aaspi_framework.h"
#include "aaspi_file.h"
#include "aaspi_segy.h"

#define MY_INT_MAX 0x7fffffff

using namespace std;

void doc()
{
	cout << "Usage aaspi_segy_writer aaspi_in= segy_out= [parameters]" << endl;
	cout << endl;
	cout << "PARAMETERS:";
	cout << "    aaspi_in" << endl;
	cout << "    Input AASPI (SEPlib/RSF) file name" << endl;
	cout << endl;
	cout << "    segy_out" << endl;
	cout << "    Output SEG-Y file name" << endl;
	cout << endl;
	cout << "    verb" << endl;
	cout << "        Verbose output [0]" << endl;
	cout << endl;
	cout << "    vblock" << endl;
	cout << "        Report progress every n traces [100]" << endl;
	cout << endl;
	cout << "    output_padded_traces   " << endl;
	cout << "        Output padded traces [n]   " << endl;
	cout << "    output_dead_traces   " << endl;
	cout << "        Output dead traces [n]   " << endl;
	cout << endl;
	cout << "    header_fn" << endl;
	cout << "        If present use text header file specified to use as the text header" << endl;
	cout << endl;
	cout << "    bfile" << endl;
	cout << "        If present use binary header file specified" << endl;
	cout << endl;
	cout << "    swap_header [0]" << endl;
	cout << "        If set to other than zero will byte swap the binary header if" << endl;
	cout << "        read from an external file with bfile= option. This is provided to" << endl;
	cout << "        support SEP/SU swapped binary header files created by Segy2sep" << endl;
	cout << endl;
	cout << "    format" << endl;
	cout << "        Specify Output Data Format [1]" << endl;
	cout << "            1=IBM floating point" << endl;
	cout << "            2=32 bit integer" << endl;
	cout << "            3=16 bit integer" << endl;
	cout << "            4=Fix point with gain (OBSOLETE!)" << endl;
	cout << "            5=IEEE floating point" << endl;
	cout << "            8=8 bit integer" << endl;
	cout << endl;
	cout << "    trmin" << endl;
	cout << "        first trace to read [1]" << endl;
	cout << endl;
	cout << "    trmax" << endl;
	cout << "        last trace to read [" <<MY_INT_MAX<<"]"<< endl;
	cout << endl;
	cout << "    ignorelist" << endl;
	cout << "        A colon separated list of names of keys to ignore for output" << endl;
	cout << endl;
	cout << "    nextra_keys" << endl;
	cout << "        Number of extra keys (Actually the number of keys to remap since" << endl;
	cout << "        we can only use the keys that we have defined in the header file." << endl;
	cout << endl;
	cout << "    extra_name#" << endl;
	cout << "        name of extra key #" << endl;
	cout << endl;
	cout << "    extra_type#" << endl;
	cout << "        data type of extra key #" << endl;
	cout << "        Types may be I8, I16, I32, IBM, IEEE, FPG" << endl;
	cout << "        I8 thru I32 are integers of the given bit length." << endl;
	cout << "        IBM & IEEE are the respective real types." << endl;
	cout << endl;
	cout << "    extra_offset#" << endl;
	cout << "        offset byte of extra key #" << endl;
	cout << "        Note: offsets are numbered starting with 1 (valid range 1-240) " << endl;
	cout << endl;
}

void fail(int xitval, const string &errmsg)
{
        cerr << "ERROR: " << errmsg << endl;
        exit(xitval);
}

list<string> process_list(const string &lstr)
{
	list<string> lst;
	int p=0;
	while (p <=lstr.size()) {
		int np=lstr.find_first_of(':',p);
		if (np == string::npos) np=lstr.size();
		string str=lstr.substr(p,np-p);
		if (! str.empty()) lst.push_back(str);
		p=np+1;
cerr << "str = " << str << " np = " << np << endl;
	}
	return lst;
}

list<string> process_extra_keys(AASPI_Framework &afw, HeaderMap &hm, int verbose)
{
	list<string> eclist;
	int nextra=afw.get_param_int("nextra_keys",0);
	if (nextra <= 0) return eclist;

	if (verbose) {
		cerr << "Remapping " << nextra << "keys" << endl;
	}
	int i;
	for (i=1; i <= nextra; ++i) {
		ostringstream oss;
		oss << "extra_name"<<i;
		if (! afw.has_param(oss.str()) ) {
			cerr << "Extra key " << i << " missing name. Skipping"
			     << endl;
			continue;
		}
		string name=afw.get_param_string(oss.str());
		                oss.str("");
		oss << "extra_offset"<<i;
		if (! afw.has_param(oss.str()) ) {
			cerr << "Extra key " << i << " missing offset. Skipping"
			     << endl;
			continue;
		}
		int offset=afw.get_param_int(oss.str(),0);
		oss.str("");
		oss << "extra_type"<<i;
		if (! afw.has_param(oss.str()) ) {
			cerr << "Extra key " << i << " missing type. Skipping"
			     << endl;
			continue;
		}
		string data_type=afw.get_param_string(oss.str());
		int dtype,dlen;
		HeaderMap::decodeDataType(data_type,dtype,dlen);
		if ( dtype == SEGY_TYPE_NONE) {
			cerr << "Invalid data type for extra key " << i
			     << " Skipping" << endl;
			continue;
		}
		if ( (offset < 1) || (offset+dlen-1 > 240) ) {
			cerr << "Invalid offset for extra key " << i
			     << " Skipping" << endl;
			continue;
		}
		if (verbose) {
			cerr << "Remapping: "<< name << " Offset="
			     << offset << " Type=" << dtype << endl;
		}
		// Sucessful remapping
		eclist.push_back(name);
		hm.add(name,offset,dtype);
	}
	return eclist;
}



void read_text_header(segyFile &sgyf, const string &th_file, int verbose)
{
	if (verbose) cerr << "Reading text header" << endl;
	segyTextHeader thead;
	ifstream ifs(th_file.c_str());
	if (ifs.fail()) fail(5,"Unable to open text header file");
	string line;
	int lcnt=0;
	while (ifs.good()) {
		getline(ifs,line);
		lcnt++;
	if (verbose) cerr << lcnt << line << endl;
		if (lcnt > 40) break;
		thead.putLine(lcnt,line);
	}
	ifs.close();
	sgyf.textHeader()=thead;
}

void build_text_header(segyFile &sgyf, AASPI_File *inf, int verbose)
{
	if (verbose) cerr << "Building text header" << endl;
	segyTextHeader thead;
	int i;
	for (i=1; i<=40; ++i) {	
		ostringstream oss;
		oss << "thead_c";
		oss.width(2); oss.fill('0');
		oss << i;
		if (inf->hasParam(oss.str())) {
			string line;
			inf->getParam(oss.str(),line);
			thead.putLine(i,line);
		} else {
			oss.str("");
			oss << "C";
			oss.width(2); oss.fill('0');
			oss << i;
			thead.putLine(i,oss.str());
		}
	}
	sgyf.textHeader()=thead;
}



void read_bin_header(segyFile &sgyf, const string &bh_file, int format, int ns, 
                     int swap_header, int verbose)
{
	if (verbose) cerr << "Reading binary header" << endl;
	// Need to add checks to see if that we create a valid SEGY file
	// Somewhere
	segyBinaryHeader bhead;
	ifstream ifs(bh_file.c_str());
	ifs.read(bhead.rawHeader(),SEGY_BIN_HEADER_SIZE);
	if (ifs.fail()) fail(6,"Unable to open binary header file");
	ifs.close();
	if (swap_header) {
		segyBinaryHeader bhead_bs;
		HeaderMap bhm;
		bhm.setBinaryHeaderDefault();
		list<string> hlist=bhm.getNameList();
		list<string>::iterator it;
		for (it=hlist.begin(); it != hlist.end(); ++it) {
			int val=bhead.getValue(*it);
			bhead_bs.putValue(*it,val);
		}
	} else {
		sgyf.binaryHeader()=bhead;
	}
}

void build_bin_header(segyFile &sgyf, AASPI_File *inf, int format, int ns, int verbose)
{
	if (verbose) cerr << "Building binary header" << endl;
	segyBinaryHeader bhead;
	int i;
	for (i=0; i<bhead.numHeaders(); ++i) {
		string nm=bhead.keyName(i);
		ostringstream oss;
		oss << "bhead_" << nm;
		int val;
		if (inf->hasParam(oss.str())) {
			inf->getParam(oss.str(),val);
			bhead.putValue(nm,val);
			if (verbose) cerr << "   " << nm << "=" << val << endl;
		}
	}
	sgyf.binaryHeader()=bhead;
}

void make_headers_consistent(segyFile &sgyf, AASPI_File *inf,
                             int ns, int format, float dt,
                             int verbose)
{
        int    output_sample_interval;
	string unit1,unit2;
	inf->getParam("unit1",unit1);
	inf->getParam("unit2",unit2);
        cerr << "unit2 = " << unit2 << endl;
        cerr << "unit1 = " << unit1 << endl;
        if (unit1 == "s" || unit1 == "kft" || unit1 == "km") {
//________________________________________________________________________________________
//         output sample is in microseconds, microkft, or microkm
//________________________________________________________________________________________
	   output_sample_interval=int(dt*1.e6+0.5);
        } else if (unit1 == "ms" || unit1 == "m" || unit1 == "ft") {
//________________________________________________________________________________________
//         output sample is in m or ft         
//________________________________________________________________________________________
	   output_sample_interval=int(dt+.5);
        } else {
           cerr << "Error! output unit1 must be s, km, ft, ms, m, or ft" << endl;
           output_sample_interval=0;  
        }

	int val;
	segyBinaryHeader bhead=sgyf.binaryHeader();
	// if job id  is zero set to a default of 1
	val=bhead.getValue("jobid");
	if (val == 0) bhead.putValue("jobid",1);
	// if line number  is zero set to a default of 1
	val=bhead.getValue("lino");
	if (val == 0) bhead.putValue("lino",1);
	// if reel number  is zero set to a default of 1
	val=bhead.getValue("reno");
	if (val == 0) bhead.putValue("reno",1);
	// if number of traces per ensemble is zero set to a default of 1
	val=bhead.getValue("ntrpr");
	if (val == 0) bhead.putValue("ntrpr",1);
	// Make sure that the sample interval is correct
	val=bhead.getValue("hdt");
	if (output_sample_interval != val) {
		cerr << "Correcting header: output_sample interval = " << output_sample_interval << unit1 << endl;
		bhead.putValue("hdt",output_sample_interval);
	}
	// if original sample interval is 0 set it to be the same as the sample interval
	val=bhead.getValue("dto");
	if (val == 0) bhead.putValue("dto",output_sample_interval);
	// Make sure that the number of samples per record is okay
	val=bhead.getValue("hns");
	if (ns != val) {
		cerr << "Correcting header: number of samples = " << ns << endl;
		bhead.putValue("hns",ns);
	}
	// if original number of samples is 0 set it to be the same as the number of samples
	val=bhead.getValue("nso");
	if (val == 0) bhead.putValue("nso",ns);
	// Make sure that the data format is correct
	val=bhead.getValue("format");
	if (format != val) {
		cerr << "Correcting header: data format = " << format << endl;
		bhead.putValue("format",format);
	}
	// If fold is 0 set to default of 1
	val=bhead.getValue("fold");
	if (val == 0) bhead.putValue("fold",1);
	// If trace sorting code is 0 set to default of 4 (stacked data)
	val=bhead.getValue("tsort");
	if (val == 0) bhead.putValue("tsort",4);
	// If vertical sum code is 0 set to default 1
	val=bhead.getValue("vscode");
	if (val == 0) bhead.putValue("vscode",1);
	// If amplitude recovery method is zero set to default 1 (none)
	val=bhead.getValue("rcvm");
	if (val == 0) bhead.putValue("rcvm",1);
	// Set measurement system code.
	if (unit2 == "m") bhead.putValue("mfeet",1);
	if (unit2 == "ft") bhead.putValue("mfeet",2);
	
	sgyf.binaryHeader()=bhead;
}

// returns true if it passes the test
bool check_segy_mapping(segyFile &sgyf, list<string> &hlist, list<string> &ek_list)
{
cerr << "Checking SEGY Mapping" << endl;
	int rv=true;
	int mapped [SEGY_TRC_HEADER_SIZE];
	int j;
	for (j=0; j<SEGY_TRC_HEADER_SIZE; ++j) mapped[j]=0;
	HeaderMap hm=sgyf.traceInterpreter().headerMap();
	vector<string> names(hlist.size());
	list<string>::iterator it;
	list<string>  rm_key_list;
	int current=1;
	for (it=hlist.begin(); it != hlist.end(); ++it) {
		header_map_entry hme=hm.get(*it);
		names[current-1]=hme.hdr_name;
		int dlen=HeaderMap::getDataTypeLength(hme.data_type);
		bool first=true;
		for (j=0; j<dlen; ++j) {
			int p=hme.start_pos-1+j;
			if (mapped[p] && first) {
				cerr << "Warning: SEGY output collision at byte "
				     << (p+1) << " between header "
				     << hme.hdr_name << " and "
				     << names[mapped[p]-1] << endl;
				list<string>::iterator fit=find(ek_list.begin(),ek_list.end(),hme.hdr_name);
				if (fit != ek_list.end()) {
					cerr << "Since key "<< hme.hdr_name << " is remapped. Removing " << names[mapped[p]-1] << endl;
					rm_key_list.push_back(names[mapped[p]-1]);
				} else {
					rv=false;
				}
				first=false;
			} else {
				mapped[p]=current;
			}
		}
		current++;
	}


	// remap keys in remove list to nowhere
	for (it=rm_key_list.begin(); it != rm_key_list.end(); ++it) {
		hm.add(*it,0,SEGY_TYPE_NONE);
	}
        // and update header map
	sgyf.traceInterpreter().setHeaderMap(hm);


	return rv;
}


void setup_trace_headers(AASPI_Framework &afw, AASPI_File *inf, segyFile& sgyf,
                         bool allow_ns_and_dt,int verbose, int& trid_key)
{

	// Get keys to ignore
	list<string> ignore;
	set<string> ignore_set;
	if (! allow_ns_and_dt) {
		ignore_set.insert("ns");
		ignore_set.insert("dt");
	}
	if (afw.has_param("ignorelist")) {
		ignore=process_list(afw.get_param_string("ignorelist"));
		list<string>::iterator it;
		for (it=ignore.begin(); it != ignore.end(); ++it) {
			ignore_set.insert(*it);
		}
	}

	// Get header map modifications from the input file
	HeaderMap hm=sgyf.traceInterpreter().headerMap();
	int sgy_nkeys=0;
	inf->getParam("segy_num_keys",sgy_nkeys);
	if (verbose) cerr << "Number of SEG-Y headers: " << sgy_nkeys << endl;

	int i;
	for (i=1; i<=sgy_nkeys; ++i) {
		ostringstream oss;
		oss << "segy_key_"<<i;
		string entry;
		inf->getParam(oss.str(),entry);
		string key;
		int off,dtyp;
		int p=0;
		int p2=entry.find(',');
		key=entry.substr(p,p2);
		p2++;
		p=entry.find(',',p2);
		istringstream iss(entry.substr(p2,p-p2));
		iss>>off;
		iss.str(entry.substr(p+1)); iss.clear();
		iss>>dtyp;

		if (ignore_set.count(key)) {  // If key is ignored set its mapping to none
			off=0;
			dtyp=SEGY_TYPE_NONE;
		}
		
		cerr << "i="<<i<< " key="<<key<<" off="<<off<<" dtyp="<<dtyp<<endl;
		hm.add(key,off,dtyp);
	}

	// Process extra keys -- not really extra keys here
	// This really allows keys to be remapped in the SEGY file from the ones
	// specified in the header file.  
	list<string> extra_key_list=process_extra_keys(afw,hm,verbose);

	sgyf.traceInterpreter().setHeaderMap(hm);

	// Set up headers to be written to the file
	list<string> klist;
	int nkey;
	inf->get_number_keys(nkey);
	for (i=1; i<=nkey; i++) {
		string name;
		inf->get_key_name(i,name);
		klist.push_back(name);
		if ( name == "trid" ) {
		     trid_key=i;
		     cerr << "trid_key = " << trid_key << endl;
		}
	}
	if ( trid_key == 0) {
	   cerr << "error! no trid key found in header format *.H@@ file!" << endl;
	   exit(1);
	}

	if (!check_segy_mapping(sgyf,klist,extra_key_list)) {
		fail(8,"Fix SEG-Y trace header mapping");
	}

	sgyf.traceInterpreter().setHeaderList(klist);
}




int main(int argc, char **argv)
{
	AASPI_Framework afw;
	if (argc < 2) {
		doc();
		exit(1);
	}
	afw.set_argvec(argc,argv);
	afw.init();

	string in_file=afw.get_param_string("aaspi_in");
	string out_file=afw.get_param_string("segy_out");
	int verbose=afw.get_param_int("verb",0);
	int vblock=afw.get_param_int("vblock",100);
	string th_file=afw.get_param_string("header_fn");
	string bh_file=afw.get_param_string("bfile");
	int swap_header=afw.get_param_int("swap_header",0);
	int format=afw.get_param_int("format",1); // Output format
	int trmin=afw.get_param_int("trmin",1);
	int trmax=afw.get_param_int("trmax",MY_INT_MAX);
	bool output_padded_traces=afw.get_param_boolean("output_padded_traces",false);
	bool output_dead_traces=afw.get_param_boolean("output_dead_traces",false);

	if (in_file.empty()) fail(1,"No input file given");
        if (out_file.empty()) fail(2,"No SEG-Y output file given");
	cout << "    header_fn = " << th_file << endl;
	
	AASPI_File *inf=afw.get_file_in(in_file);
	if (!inf) fail(3,"Uable to open input file");

	int trid_key;
	int nax=inf->get_num_axes();
	int ns;
	float t0,dt;
	string lab;
	inf->get_axis_par(1,ns,t0,dt,lab);
	if (verbose) {
		cerr << "Number of axes:   " << nax << endl;
		cerr << "Samples/trace:    " << ns << endl;
		cerr << "Sample interval:  " << dt*1000.f << " ms" << endl;
		cerr << "Start time:       " << t0 << " s" << endl;
		cerr << "End time:         " << (t0+(ns-1)*dt) << " s" << endl;
	}
        cerr << " output padded traces = " << output_padded_traces << endl;
        cerr << " output dead   traces = " << output_dead_traces << endl;
	int i,ntraces=1;
	for (i=2; i<=nax; ++i) {
		int n;
		float o,d;
		string label;
		inf->get_axis_par(i,n,o,d,label);
		ntraces *= n;
		if (verbose) {
			cerr << "  Axis " << i << " n=";
			cerr.width(10);
			cerr << n << "    From " ;
			cerr.width(8);
			cerr << o << " to ";
			cerr.width(8);
			cerr << (o+(n-1)*d) << " by ";
			cerr.width(8);	
			cerr << d << "   " << label << endl;
		}
	}
	int esize=4;
        string unit1;
        string unit2;
	inf->getParam("esize",esize);
	inf->getParam("unit1",unit1);
	inf->getParam("unit2",unit2);
        cerr << " after getParam. unit1 = " << unit1 << endl;
        cerr << " after getParam. unit2 = " << unit2 << endl;
	if (verbose) {
		cerr << "Total traces: " << ntraces << endl;
	}

	// Indicate and fix data format
	if (format != 1) {
		cerr << "Data format requested=" << format;
		switch (format) {
			case SEGY_TYPE_IEEE:
				cerr << " (IEEE float)" << endl;
				break;
			case SEGY_TYPE_I32:
				cerr << " (32 bit integer)" << endl;
				break;
			case SEGY_TYPE_I16:
				cerr << " (16 bit integer)" << endl;
				break;
			case SEGY_TYPE_I8:
				cerr << " (8 bit integer)" << endl;
				break;
			case SEGY_TYPE_FPG:
				cerr << " (Fixed point with gain)" << endl;
				break;
			default:
				cerr << " (UNKNOWN) Setting format=1 (IBM)" << endl;
				format=1;
		}
	}
	if (verbose) {
		cerr << "Using data format =" << format;
		switch (format) {
			case SEGY_TYPE_IBM:
				cerr << " (IBM float)" << endl;
				break;
			case SEGY_TYPE_IEEE:
				cerr << " (IEEE float)" << endl;
				break;
			case SEGY_TYPE_I32:
				cerr << " (32 bit integer)" << endl;
				break;
			case SEGY_TYPE_I16:
				cerr << " (16 bit integer)" << endl;
				break;
			case SEGY_TYPE_I8:
				cerr << " (8 bit integer)" << endl;
				break;
			case SEGY_TYPE_FPG:
				cerr << " (Fixed point with gain)" << endl;
		}
	}
	

	segyFile sgyf;
	bool rv=sgyf.OpenOut(out_file,ns,dt,format,true);
	if (!rv) fail(4,"Uable to open output file");

	if (afw.has_param("header_fn")) {
		read_text_header(sgyf,th_file,verbose);
	} else {
		build_text_header(sgyf,inf,verbose);
	}

	if (afw.has_param("bfile")) {
		read_bin_header(sgyf,bh_file,format,ns,swap_header,verbose);
	} else {
		build_bin_header(sgyf,inf,format,ns,verbose);
	}
	make_headers_consistent(sgyf,inf,ns,format,dt,verbose);
	
	rv=sgyf.write_headers();
	if (!rv) fail(7,"Error writing SEG-Y file headers");

	// Set up trace_header mappings
	bool allow_ns_and_dt=false;
	setup_trace_headers(afw,inf,sgyf,allow_ns_and_dt,verbose,trid_key);
//      cerr << "after setup_trace_headers trid_key = " << trid_key;

	// Set laga override
	int laga;
//	string unit1,unit2;
	inf->getParam("unit1",unit1);
	inf->getParam("unit2",unit2);
        if (unit1 == "s" || unit1 == "kft" || unit1 == "km") {
	   if (t0 >= 0.f) { // nint
		laga=int(1000.0f*t0+0.5f);
	   } else {
		laga=-int(-1000.0f*t0+0.5f);
	   }
	   if (verbose) cerr << "Setting lagA to " << laga << " milli"+unit1 << endl;
        } else if (unit1 == "ms" || unit1 == "ft" || unit1 == "m") {
	   if (t0 >= 0.f) { // nint
		laga=int(1.0f*t0+0.5f);
	   } else {
		laga=-int(-1.0f*t0+0.5f);
	   }
	   if (verbose) cerr << "Setting lagA to " << laga <<  unit1 << endl;
        } else {
           cerr << "Error! output unit1 must be s, km, ft, ms, m, or ft" << endl;
           laga=0;
        }
	// lagA (16 bit int at byte 105)
	sgyf.traceInterpreter().addOverride(105,SEGY_TYPE_I16,laga);

	// Thang Ha, 20170113: Setting delrt to be laga.
	// Needed for Petrel to recognize the starting time.
	if (verbose) {
		cerr << "Setting delrt to " << laga  << endl;
	}
	// delrt (16 bit int at byte 109)
	sgyf.traceInterpreter().addOverride(109, SEGY_TYPE_I16, laga);

	int nkeys;
	inf->get_number_keys(nkeys);
	int *kbuff=new int [nkeys];

	int rsize=esize*ns;
	float *dbuff= new float [ns];

	// Now write the traces to the file

	if (trmax > ntraces) trmax=ntraces;
	aaspi_int64 pos=esize*ns*aaspi_int64(trmin-1);
	inf->Seek(pos);
	int trn, tc=0;
	for (trn=trmin; trn<=trmax; ++trn) {
		inf->get_val_headers(trn,1,kbuff);
 //		int i;
 // 		for (i=0; i<nkeys; ++i) {
 // 			cerr << i<< " " <<kbuff[i]<< endl;
 // 		}
		inf->Read(dbuff,rsize);
//                 cerr << " trid_key ="<< trid_key <<" kbuff[trid_key-1] " << kbuff[trid_key-1] <<  endl;
                if ( kbuff[trid_key-1] == 1 ) {
//                 cerr << " write out the live trace"  <<  endl;
//_____________________________________________________________
//                  a live trace. Always write these out!
//_____________________________________________________________
   	           sgyf.writeTraces(1,ns,nkeys,dbuff,kbuff);
                } else if ( kbuff[trid_key-1] == 2 && output_dead_traces ) {
//_____________________________________________________________
//                 we want to write out dead traces as well
//_____________________________________________________________
   	           sgyf.writeTraces(1,ns,nkeys,dbuff,kbuff);
//                 cerr << " trid_key ="<< trid_key <<" kbuff[trid_key-1] " << kbuff[trid_key-1] << "  write dead trace" << endl;
                } else if ( kbuff[trid_key-1] == 3 && output_padded_traces) {
//_____________________________________________________________
//                 we also wish to write out padded traces
//_____________________________________________________________
   	           sgyf.writeTraces(1,ns,nkeys,dbuff,kbuff);
                }
		tc++;
		if (tc%vblock == 0) cerr << tc << " traces written" << endl;
	}
	cerr << tc << " traces written" << endl;

	delete [] dbuff;

	if (verbose) cerr << "Closing output SEG-Y file" << endl;
	sgyf.Close();
	if (verbose) cerr << "Closing input file" << endl;
	inf->Close();

	return 0;
}