handling PRF change for ALOS-1 raw data

LT1AB
CunrenLiang 2020-10-19 19:45:01 -07:00 committed by CunrenLiang
parent e9bd7edeb3
commit 9837bf381c
22 changed files with 2390 additions and 115 deletions

View File

@ -25,7 +25,12 @@
# Author: Walter Szeliga
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
#********************************************************************************
#This program has been upgraded to handle the ALOS-1 PRF change issue.
#
#Cunren Liang, 12-December-2019
#California Institute of Technology, Pasadena, CA
#********************************************************************************
import os
@ -114,8 +119,8 @@ class ALOS(Sensor):
complex(-6.297074e-3,8.026685e-3),
complex(7.217117e-1,-2.367683e-2))
constants = Constants(iBias=15.5,
qBias=15.5,
constants = Constants(iBias=63.5,
qBias=63.5,
pointingDirection=-1,
antennaLength=8.9)
@ -499,7 +504,7 @@ class ALOS(Sensor):
outputNow = self.output + appendStr
if not (self._resampleFlag == ''):
filein = self.output + '__tmp__'
self.imageFile.extractImage(filein)
self.imageFile.extractImage(filein, i) #image number start with 0
self.populateMetadata()
objResample = None
if(self._resampleFlag == 'single2dual'):
@ -513,7 +518,7 @@ class ALOS(Sensor):
objResample.updateFrame(self.frame)
os.remove(filein)
else:
self.imageFile.extractImage(outputNow)
self.imageFile.extractImage(outputNow, i) #image number start with 0
self.populateMetadata()
width = self.frame.getImage().getWidth()
# self.readOrbitPulse(self._leaderFile,outputNow,width)
@ -721,7 +726,7 @@ class ImageFile(object):
return None
def extractImage(self,output=None):
def extractImage(self,output=None, image_i=0):
"""For now, call a wrapped version of ALOS_pre_process"""
productLevel = float(self.parent.leaderFile.sceneHeaderRecord.metadata[
'Product level code'])
@ -731,13 +736,13 @@ class ImageFile(object):
elif productLevel == 1.1:
self.extractSLC(output)
elif productLevel == 1.0:
self.extractRaw(output)
self.extractRaw(output, image_i) #image number start with 0
else:
raise ValueError(productLevel)
return None
@use_api
def extractRaw(self,output=None):
def extractRaw(self,output=None, image_i=0):
#if (self.numberOfSarChannels == 1):
# print "Single Pol Data Found"
# self.extractSinglePolImage(output=output)
@ -748,15 +753,16 @@ class ImageFile(object):
if self.parent.leaderFile.sceneHeaderRecord.metadata[
'Processing facility identifier'] == 'ERSDAC':
prmDict = alos.alose_Py(self.parent._leaderFile,
self.parent._imageFile, output)
self.parent._imageFile, output, image_i) #image number start with 0
else:
prmDict = alos.alos_Py(self.parent._leaderFile,
self.parent._imageFile, output)
self.parent._imageFile, output, image_i) #image number start with 0
pass
# updated 07/24/2012
self.width = prmDict['NUMBER_BYTES_PER_LINE'] - 2 * prmDict['FIRST_SAMPLE']
self.length = self.imageFDR.metadata['Number of lines per data set']
#self.length = self.imageFDR.metadata['Number of lines per data set']
self.length = prmDict['NUMBER_LINES']
self.prefix = self.imageFDR.metadata[
'Number of bytes of prefix data per record']
self.suffix = self.imageFDR.metadata[

View File

@ -25,7 +25,8 @@
// Author: Giangi Sacco
//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
//program updated to handle PRF change issue of ALOS-1
//Cunren Liang, December 2019
#include <Python.h>
@ -68,11 +69,12 @@ PyInit_alos()
PyObject *alos_C(PyObject* self,PyObject *args)
{
char *imageFile,*leaderFile,*outFile;
int image_i;
struct PRM inputPRM;
struct PRM outputPRM;
struct GLOBALS globals;
if(!PyArg_ParseTuple(args,"sss",&leaderFile,&imageFile,&outFile))
if(!PyArg_ParseTuple(args,"sssi",&leaderFile,&imageFile,&outFile, &image_i))
{
return NULL;
}
@ -96,7 +98,7 @@ PyObject *alos_C(PyObject* self,PyObject *args)
globals.dopp = 0; // Are we calculating a doppler?
globals.tbias = 0.0; // Is there a time bias to fix poor orbits?
ALOS_pre_process(inputPRM,&outputPRM,globals);
ALOS_pre_process(inputPRM,&outputPRM,globals,image_i);
PyObject * dict = PyDict_New();
createDictionaryOutput(&outputPRM,dict);
@ -106,11 +108,12 @@ PyObject *alos_C(PyObject* self,PyObject *args)
PyObject *alose_C(PyObject* self,PyObject *args)
{
char *imageFile,*leaderFile,*outFile;
int image_i;
struct PRM inputPRM;
struct PRM outputPRM;
struct GLOBALS globals;
if(!PyArg_ParseTuple(args,"sss",&leaderFile,&imageFile,&outFile))
if(!PyArg_ParseTuple(args,"sssi",&leaderFile,&imageFile,&outFile, &image_i))
{
return NULL;
}
@ -134,7 +137,7 @@ PyObject *alose_C(PyObject* self,PyObject *args)
globals.dopp = 0; // Are we calculating a doppler?
globals.tbias = 0.0; // Is there a time bias to fix poor orbits?
ALOS_pre_process(inputPRM,&outputPRM,globals);
ALOS_pre_process(inputPRM,&outputPRM,globals,image_i);
PyObject * dict = PyDict_New();
createDictionaryOutput(&outputPRM,dict);
@ -184,6 +187,14 @@ PyObject * createDictionaryOutput(struct PRM * prm, PyObject * dict)
intVal = PyLong_FromLong((long) prm->good_bytes);
PyDict_SetItemString(dict,"NUMBER_GOOD_BYTES",intVal);
Py_XDECREF(intVal);
intVal = PyLong_FromLong((long) prm->num_lines);
PyDict_SetItemString(dict,"NUMBER_LINES",intVal);
Py_XDECREF(intVal);
intVal = PyLong_FromLong((long) prm->num_rng_bins);
PyDict_SetItemString(dict,"NUMBER_RANGE_BIN",intVal);
Py_XDECREF(intVal);

View File

@ -25,7 +25,8 @@
// Author: Giangi Sacco
//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
//program updated to handle PRF change issue of ALOS-1
//Cunren Liang, December 2019
#ifndef alosmodule_h
@ -42,7 +43,7 @@ extern "C"
PyObject *alose_C(PyObject *self,PyObject *args);
PyObject *createDictionaryOutput(struct PRM *prm,PyObject *dict);
int ALOS_pre_process(struct PRM inputPRM, struct PRM *outputPRM,
struct GLOBALS globals);
struct GLOBALS globals, int image_i);
}
static PyMethodDef alos_methods[] =

View File

@ -17,9 +17,25 @@
* added write_roi
* *****************************************************************************/
/********************************************************************************
This program has been upgraded to handle the ALOS-1 PRF change issue.
Cunren Liang, 12-December-2019
California Institute of Technology, Pasadena, CA
*********************************************************************************/
#include "alosglobals.h"
#include"image_sio.h"
#include"lib_functions.h"
#include "resamp.h"
//ALOS I or Q mean = 15.5, so get 15 or 16 randomly here
//#define ZERO_VALUE (char)(15 + rand() % 2)
//I changed the dynamic range when reading data
//ALOS I or Q mean = 63.5, so get 63 or 64 randomly here
#define ZERO_VALUE (char)(63 + rand() % 2)
char *USAGE = "\n\nUsage: ALOS_pre_process imagefile LEDfile [-near near_range] [-radius RE] [-swap] [-V] [-debug] [-quiet] \n"
"\ncreates data.raw and writes out parameters (PRM format) to stdout\n"
@ -44,9 +60,8 @@ char *USAGE = "\n\nUsage: ALOS_pre_process imagefile LEDfile [-near near_rang
"-tbias tbias correct the clock bias (positive value means plus)\n"
"Example:\n"
"ALOS_pre_process IMG-HH-ALPSRP050420840-H1.0__A LED-ALPSRP050420840-H1.0__A \n";
long read_ALOS_data (FILE *, FILE *, struct PRM *, long *);
long read_ALOSE_data (FILE *, FILE *, struct PRM *, long *);
long read_ALOS_data (FILE *, FILE *, struct PRM *, long *, struct resamp_info *, int);
long read_ALOSE_data (FILE *, FILE *, struct PRM *, long *, struct resamp_info *, int);
void parse_ALOS_commands(int, char **, char *, struct PRM *);
void set_ALOS_defaults(struct PRM *);
void print_ALOS_defaults(struct PRM *);
@ -58,19 +73,54 @@ int write_roi(char *, FILE *, struct PRM, struct ALOS_ORB, char *);
// ISCE stuff
void init_from_PRM(struct PRM inPRM, struct PRM *prm);
int resamp_azimuth(char *slc2, char *rslc2, int nrg, int naz1, int naz2, double prf, double *dopcoeff, double *azcoef, int n, double beta);
int
ALOS_pre_process(struct PRM inputPRM, struct PRM *outputPRM,struct GLOBALS globals)
ALOS_pre_process(struct PRM inputPRM, struct PRM *outputPRM,struct GLOBALS globals, int image_i) //image number starts with 0!!!
{
FILE *imagefile, *ldrfile;
FILE *rawfile[11];//*prmfile[11];
//char prmfilename[128];
FILE *rawfile[11], *prmfile[11];
char prmfilename[128];
int nPRF;
long byte_offset;
struct PRM prm;
struct ALOS_ORB orb;
char date[8];
//////////////////////////////////////////////
FILE *resampinfofile;
struct resamp_info rspi;
struct resamp_info rspi_new;
struct resamp_info rspi_pre[100];//maximum number of frames: 100
int i, j, k;
double SC_clock_start;
double SC_clock_start_resamp;
double d2s = 24.0 * 3600.0;
double line_number_first;
int num_lines_out;
int gap_flag;
double prf_all[200];//maximum number of prfs: 200
int frame_counter_start_all[200];//maximum number of prfs: 200
int nPRF_all;//maximum number of prfs: 200
double dopcoeff[4];
double azcoef[2];
int num_lines_max, j_max;
char outputfile[256];
char *data;
FILE *first_prf_fp;
FILE *next_prf_fp;
int num_lines_append;
//int num_lines_gap;
int ret;
//////////////////////////////////////////////
//if (argc < 3) die (USAGE,"");
printf("reading image: %d\n", image_i);
/* set flags */
dopp = globals.dopp;
@ -91,19 +141,21 @@ char date[8];
init_from_PRM(inputPRM,&prm);
//parse_ALOS_commands(argc, argv, USAGE, &prm);
//if (verbose) print_ALOS_defaults(&prm);
/* apply an additional timing bias based on corner reflector analysis */
//tbias = tbias - 0.0020835;
if (verbose) print_ALOS_defaults(&prm);
if (is_big_endian_() == -1) {swap = 1;fprintf(stderr,".... swapping bytes\n");} else {swap = 0;}
/* IMG and LED files should exist already */
if ((rawfile[0] = fopen(prm.input_file,"w")) == NULL) die("can't open ",prm.input_file);
if ((imagefile = fopen(globals.imagefilename, "r")) == NULL) die ("couldn't open Level 1.0 IMG file \n",globals.imagefilename);
if ((ldrfile = fopen(inputPRM.led_file, "r")) == NULL) die ("couldn't open LED file \n",inputPRM.led_file);
/* if it exists, copy to prm structure */
//strcpy(prm.led_file,leaderFilename);
strcpy(prm.led_file,inputPRM.led_file);
/* name and open output files and header files for raw data (but input for later processing) */
//get_files(&prm, &rawfile[nPRF], &prmfile[nPRF], prmfilename, argv[1], nPRF);
get_files(&prm, &rawfile[nPRF], &prmfile[nPRF], prmfilename, prm.input_file, nPRF);
/* read sarleader; put info into prm; write log file if specified */
read_ALOS_sarleader(ldrfile, &prm, &orb);
@ -125,7 +177,7 @@ char date[8];
/* if prf changes, create new prm and data files */
if (nPRF > 0 ) {
if (verbose) fprintf(stderr,"creating multiple files due to PRF change (*.%d) \n",nPRF+1);
//get_files(&prm, &rawfile[nPRF], &prmfile[nPRF], prmfilename, argv[1], nPRF);
get_files(&prm, &rawfile[nPRF], &prmfile[nPRF], prmfilename, prm.input_file, nPRF);
}
/* set the chirp extension to 500 if FBD fs = 16000000 */
@ -141,21 +193,26 @@ char date[8];
returns byte offset if the PRF changes */
/* calculate parameters from orbit */
if (ALOS_format == 0) {
byte_offset = read_ALOS_data(imagefile, rawfile[nPRF], &prm, &byte_offset);
byte_offset = read_ALOS_data(imagefile, rawfile[nPRF], &prm, &byte_offset, &rspi, nPRF);
}
/* ERSDAC - use read_ALOSE_data */
if (ALOS_format == 1) {
byte_offset = read_ALOSE_data(imagefile, rawfile[nPRF], &prm, &byte_offset);
byte_offset = read_ALOSE_data(imagefile, rawfile[nPRF], &prm, &byte_offset, &rspi, nPRF);
}
//fclose(rawfile[nPRF]);
// should work for AUIG and ERSDAC
ALOS_ldr_orbit(&orb, &prm);
/* calculate doppler from raw file */
dopp=1;//always compute doppler for doing prf resampling
if (dopp == 1) calc_dop(&prm);
//prf as a function of range in Hz
rspi.fd1[nPRF] = prm.fd1;
rspi.fdd1[nPRF] = prm.fdd1;
rspi.fddd1[nPRF] = prm.fddd1;
//rspi.input_file[nPRF] = prm.input_file;
strcpy(rspi.input_file[nPRF], prm.input_file);
/* divide prf in half for quad_pol */
/* fix chirp slope */
@ -172,7 +229,7 @@ char date[8];
if (force_slope == 1) prm.chirp_slope = forced_slope;
/* write ascii output, SIO format */
//put_sio_struct(prm, prmfile[nPRF]);
put_sio_struct(prm, prmfile[nPRF]);
/* write roi_pac output */
if (roi) {
@ -184,6 +241,322 @@ char date[8];
nPRF++;
}
rspi.nPRF=nPRF;
////////////////////////////////////////////////////////////////////////////////////////////////////////////////
printf("\nPRF details of frame: %d\n", image_i);
printf("+++++++++++++++++++++++++++++++++++++++++++++++\n");
printf("number of PRF: %d\n", rspi.nPRF);
for (i = 0; i < rspi.nPRF; i++){
printf("PRF %d prf (Hz): %f\n", i+1, rspi.prf[i]);
printf("PRF %d start time (days): %20.12f\n", i+1, rspi.SC_clock_start[i]);
printf("PRF %d frame_counter_start: %d\n", i+1, rspi.frame_counter_start[i]);
printf("PRF %d frame_counter_end: %d\n", i+1, rspi.frame_counter_end[i]);
printf("PRF %d number of lines: %d\n\n", i+1, rspi.frame_counter_end[i]-rspi.frame_counter_start[i]+1);
}
//open parameter file for doing time adjustment and interpolation
if (image_i == 0){
if((resampinfofile = fopen("resampinfo.bin", "wb")) == NULL)
die("couldn't open resampinfo file","resampinfo.bin");
}
else{
//open the file for reading and appending
if((resampinfofile = fopen("resampinfo.bin", "ab+")) == NULL)
die("couldn't open resampinfo file","resampinfo.bin");
rewind(resampinfofile);
for(i=0; i < image_i; i++){
if((fread((void *) &rspi_pre[i],sizeof(struct resamp_info), 1, resampinfofile)) != 1)
die("couldn't read from file","resampinfo.bin");
}
}
//get parameter from this image
memcpy(&rspi_pre[image_i], &rspi, sizeof(struct resamp_info));
//initialize rspi_new with resamp_info from reading the image, put the adjusted time in it
memcpy(&rspi_new, &rspi, sizeof(struct resamp_info));
//adjust start time
//unified PRF of the full track: first prf of first image
//start time of the full track: first line of first image
//only adjust time when the format is not ERSDAC format, becasue ERSDAC format does not have sdr.frame_counter.
printf("adjust start times\n");
if(ALOS_format == 0){
if(image_i==0){
//adjust start time of prf file i, no need to adjust for first prf
for(i = 1; i < rspi_pre[0].nPRF; i++){
//time of the line just before the first line of first prf file
SC_clock_start = rspi_pre[0].SC_clock_start[0] - (1.0/rspi_pre[0].prf[0]) / d2s;
//time of the last line of each prf file
for(j = 0; j < i; j++){
if(rspi_pre[0].num_lines[j] != rspi_pre[0].frame_counter_end[j] - rspi_pre[0].frame_counter_start[j] + 1)
fprintf(stderr, "\n\nWARNING: in image %d prf file %d, \
number of lines in file: %d is not equal to that computed from frame_counter: %d\n\n", \
0, j, rspi_pre[0].num_lines[j], rspi_pre[0].frame_counter_end[j] - rspi_pre[0].frame_counter_start[j] + 1);
SC_clock_start += (rspi_pre[0].frame_counter_end[j]-rspi_pre[0].frame_counter_start[j]+1) * (1.0/rspi_pre[0].prf[j]) / d2s;
}
//time of the first line of current prf file
SC_clock_start += (1.0/rspi_pre[0].prf[i]) / d2s;
printf("time adjustment result for image %d, prf %d:\n", image_i, i);
printf("+++++++++++++++++++++++++++++++++++++++++++++++\n");
printf("original start time: %20.12f\n", rspi_pre[0].SC_clock_start[i]);
printf("adjusted start time: %20.12f\n", SC_clock_start);
printf("original - adjusted: %f (number of PRI)\n\n", (rspi_pre[0].SC_clock_start[i]-SC_clock_start)*d2s/(1.0/rspi_pre[0].prf[i]));
//update
rspi_new.SC_clock_start[i] = SC_clock_start;
}
}
else{
//1. check to see if there is gap between images
gap_flag = 0;
for(i = 0; i < image_i; i++){
if (rspi_pre[i].frame_counter_end[rspi_pre[i].nPRF-1] - rspi_pre[i+1].frame_counter_start[0] <= -2){
fprintf(stderr, "\n\nWARNING: there are gaps between image %d and image: %d\n", i, i+1);
fprintf(stderr, "since we don't know the prf of these gap lines, we are not able to adjust starting time\n\n");
gap_flag = 1;
}
}
//2. adjust start time
if(gap_flag == 0){
//2.1 count the number of prf chunks in the full track including this image
nPRF_all = 0;
for(i = 0; i < image_i+1; i++){
for(j = 0; j < rspi_pre[i].nPRF; j++){
if((i==0) && (j==0)){
prf_all[nPRF_all] = rspi_pre[i].prf[j];
frame_counter_start_all[nPRF_all] = rspi_pre[i].frame_counter_start[j];
nPRF_all += 1;
}
else{
if((rspi_pre[i].frame_counter_start[j]>frame_counter_start_all[nPRF_all-1]) && (rspi_pre[i].prf[j]!=prf_all[nPRF_all-1])){
prf_all[nPRF_all] = rspi_pre[i].prf[j];
frame_counter_start_all[nPRF_all] = rspi_pre[i].frame_counter_start[j];
nPRF_all += 1;
}
}
}
}
printf("number of prfs including this image: %d\n", nPRF_all);
printf("list of prfs:\n");
for(i = 0; i < nPRF_all; i++){
printf("frame_counter: %d, prf: %f\n", frame_counter_start_all[i], prf_all[i]);
}
//2.2 adjust start time
for(i = 0; i < rspi_pre[image_i].nPRF; i++){
//time of the line just before the first line of first prf file
//because the unite is day, the errors caused can be 0.042529743164777756 lines, should remove the integer or year part of SC_clock_start, or
//use second as unit in the future
SC_clock_start = rspi_pre[0].SC_clock_start[0] - (1.0/rspi_pre[0].prf[0]) / d2s;
//if there is only one PRF (no prf changes across all images)
if(nPRF_all == 1){
SC_clock_start += (rspi_pre[image_i].frame_counter_start[0] - rspi_pre[0].frame_counter_start[0] + 1) * (1.0/rspi_pre[0].prf[0]) / d2s;
}
else{
//find its position among the prfs, start from the second prf
for(j = 1; j < nPRF_all; j++){
if(rspi_pre[image_i].frame_counter_start[i] < frame_counter_start_all[j]){
//time of the last line of each prf chuck
for(k = 1; k < j; k++)
SC_clock_start += (frame_counter_start_all[k]-frame_counter_start_all[k-1]) * (1.0/prf_all[k-1]) / d2s;
SC_clock_start += (rspi_pre[image_i].frame_counter_start[i] - frame_counter_start_all[j-1] + 1) * (1.0/prf_all[j-1]) / d2s;
break;
}
else if(rspi_pre[image_i].frame_counter_start[i] == frame_counter_start_all[j]){
//time of the last line of each prf chuck
for(k = 1; k < j; k++)
SC_clock_start += (frame_counter_start_all[k]-frame_counter_start_all[k-1]) * (1.0/prf_all[k-1]) / d2s;
SC_clock_start += (rspi_pre[image_i].frame_counter_start[i] - frame_counter_start_all[j-1] + 1) * (1.0/prf_all[j-1]) / d2s;
//extra pri of j-1 above, so remove it and add the pri of j
SC_clock_start += (1.0/prf_all[j]) / d2s - (1.0/prf_all[j-1]) / d2s;
break;
}
else{
if(j == nPRF_all - 1){
for(k = 1; k < j+1; k++)
SC_clock_start += (frame_counter_start_all[k]-frame_counter_start_all[k-1]) * (1.0/prf_all[k-1]) / d2s;
SC_clock_start += (rspi_pre[image_i].frame_counter_start[i] - frame_counter_start_all[j] + 1) * (1.0/prf_all[j]) / d2s;
break;
}
else{
continue;
}
}
}
}
//time of the first line of current prf file
printf("time adjustment result for image %d, prf %d:\n", image_i, i);
printf("+++++++++++++++++++++++++++++++++++++++++++++++\n");
printf("original start time: %20.12f\n", rspi_pre[image_i].SC_clock_start[i]);
printf("adjusted start time: %20.12f\n", SC_clock_start);
printf("original - adjusted: %f (number of PRI)\n\n", (rspi_pre[image_i].SC_clock_start[i]-SC_clock_start)*d2s/(1.0/rspi_pre[image_i].prf[i]));
//update
rspi_new.SC_clock_start[i] = SC_clock_start;
}
}
}
}
// use parameters from rspi_pre[image_i], instead of rspi_new (to be updated)
//except rspi_new.SC_clock_start[i], since it was updated (more accurate) above.
printf("azimuth resampling\n");
for(i = 0; i < rspi_pre[image_i].nPRF; i++){
if((image_i==0)&&(i==0))
continue;
//convention: line numbers start with zero
//line number of first line of first prf of first image: 0
//line number of first line of this prf file
line_number_first = (rspi_new.SC_clock_start[i] - rspi_pre[0].SC_clock_start[0]) * d2s / (1.0 / rspi_pre[0].prf[0]);
//unit: pri of first prf of first image
num_lines_out = (rspi_pre[image_i].frame_counter_end[i] - rspi_pre[image_i].frame_counter_start[i] + 1) * (1.0/rspi_pre[image_i].prf[i]) / (1.0/rspi_pre[0].prf[0]);
if((fabs(roundfi(line_number_first)-line_number_first)<0.1) && (rspi_pre[image_i].prf[i]==rspi_pre[0].prf[0]))
continue;
//time of first line of the resampled image
SC_clock_start_resamp = rspi_pre[0].SC_clock_start[0] + roundfi(line_number_first) * (1.0 / rspi_pre[0].prf[0]) / d2s;
//compute offset parameters
//azcoef[0] + azpos * azcoef[1]
azcoef[0] = (SC_clock_start_resamp - rspi_new.SC_clock_start[i]) * d2s / (1.0/rspi_pre[image_i].prf[i]);
azcoef[1] = (1.0/rspi_pre[0].prf[0]) / (1.0/rspi_pre[image_i].prf[i]) - 1.0;
//use doppler centroid frequency estimated from prf with maximum number of lines in this image
num_lines_max = -1;
j_max = -1;
for(j = 0; j < rspi_pre[image_i].nPRF; j++){
if(rspi_pre[image_i].num_lines[j] >= num_lines_max){
num_lines_max = rspi_pre[image_i].num_lines[j];
j_max = j;
}
}
dopcoeff[0] = rspi_pre[image_i].fd1[j_max]; //average prf for alos-1 is good enough (calc_dop.c).
dopcoeff[1] = 0.0;
dopcoeff[2] = 0.0;
dopcoeff[3] = 0.0;
//The filenames of all three files created for each prf, are from prm.input_file
//PRM: prm.input_file.PRM + (.prfno_start_from_1, if not first prf)
//data: prm.input_file + (.prfno_start_from_1, if not first prf)
//data after resampling: prm.input_file + (.prfno_start_from_1, if not first prf) + .interp
sprintf(outputfile,"%s.interp", rspi_pre[image_i].input_file[i]);
//start interpolation
resamp_azimuth(rspi_pre[image_i].input_file[i], outputfile, rspi_pre[image_i].num_bins[i], num_lines_out, rspi_pre[image_i].num_lines[i], rspi_pre[image_i].prf[i], dopcoeff, azcoef, 9, 5.0);
//update parameters
rspi_new.SC_clock_start[i] = SC_clock_start_resamp;
rspi_new.num_lines[i] = num_lines_out;
rspi_new.prf[i] = rspi_pre[0].prf[0];
rspi_new.fd1[i] = dopcoeff[0];
rspi_new.fdd1[i]= dopcoeff[1];
rspi_new.fddd1[i]=dopcoeff[2];
strcpy(rspi_new.input_file[i], outputfile);
}
//concatenate prfs: put all prfs to the first prf
// use parameters from rspi_new (updated), instead of rspi_pre[image_i]
if(rspi_new.nPRF > 1){
//prepare for appending subsequent prfs to first prf: open files and allocate memory
if((first_prf_fp = fopen(rspi_new.input_file[0], "ab")) == NULL)
die("can't open", rspi_new.input_file[0]);
//number of range samples in each prf is asummed to be same
if((data = (char *)malloc(2*sizeof(char)*rspi_new.num_bins[0])) == NULL)
die("can't allocate memory for data", "");
//append prf i
for(i = 1; i < rspi_new.nPRF; i++){
//number of lines to be appended between frames if there are gaps
num_lines_append = (rspi_new.SC_clock_start[i] - rspi_new.SC_clock_start[0]) * d2s / (1.0/rspi_pre[0].prf[0]) - rspi_new.num_lines[0];
if(num_lines_append >= 1){
for(j = 0; j < num_lines_append; j++){
for(k = 0; k < 2*rspi_new.num_bins[i]; k++)
data[k] = ZERO_VALUE;
if(fwrite((char *)data, 2*sizeof(char)*rspi_new.num_bins[i], 1, first_prf_fp) != 1)
die("can't write data to", rspi_new.input_file[0]);
}
rspi_new.num_lines[0] += num_lines_append;
}
//append data from rspi_new.input_file[i]
if((next_prf_fp = fopen(rspi_new.input_file[i], "rb")) == NULL)
die("can't open", rspi_new.input_file[i]);
num_lines_append = 0;
for(j = 0; j < rspi_new.num_lines[i]; j++){
if((rspi_new.SC_clock_start[i] + j * (1.0/rspi_pre[0].prf[0]) / d2s - rspi_new.SC_clock_start[0]) * d2s / (1.0/rspi_pre[0].prf[0]) >= rspi_new.num_lines[0]){
if(fread((char *)data, 2*sizeof(char)*rspi_new.num_bins[i], 1, next_prf_fp) != 1)
die("can't read data from", rspi_new.input_file[i]);
if(fwrite((char *)data, 2*sizeof(char)*rspi_new.num_bins[i], 1, first_prf_fp) != 1)
die("can't write data to", rspi_new.input_file[0]);
num_lines_append += 1;
}
else{
fseek(next_prf_fp, 2*sizeof(char)*rspi_new.num_bins[i], SEEK_CUR);
}
}
rspi_new.num_lines[0] += num_lines_append;
fclose(next_prf_fp);
}
free(data);
fclose(first_prf_fp);
}
//tidy up intermediate files
for(i = 0; i < rspi_pre[image_i].nPRF; i++){
//if Return value = 0 then it indicates str1 is equal to str2.
ret = strcmp(rspi_new.input_file[i], rspi_pre[image_i].input_file[i]);
if(i == 0){
if(ret != 0){
//remove original
if(remove(rspi_pre[image_i].input_file[i]) != 0)
die("can't delete file", rspi_pre[image_i].input_file[i]);
//keep resampled and appended
if(rename(rspi_new.input_file[i], rspi_pre[image_i].input_file[i]) != 0)
die("can't rename file", rspi_new.input_file[i]);
}
}
else{
//remove original
if(remove(rspi_pre[image_i].input_file[i]) != 0)
die("can't delete file", rspi_pre[image_i].input_file[i]);
//remove resampled
if(ret != 0){
if(remove(rspi_new.input_file[i]) != 0)
die("can't delete file", rspi_new.input_file[i]);
}
}
}
//update prm
prm.prf = rspi_new.prf[0];
prm.num_lines = rspi_new.num_lines[0];
prm.SC_clock_start = rspi_new.SC_clock_start[0];
prm.SC_clock_stop = prm.SC_clock_start + (prm.num_lines - 1) * (1.0/prm.prf) / d2s;
prm.fd1 = rspi_pre[image_i].fd1[j_max]; //average prf for alos-1 is good enough (calc_dop.c).
prm.fdd1 = 0.0;
prm.fddd1 =0.0;
prm.xmi = 63.5;
prm.xmq = 63.5;
//write to resampinfo.bin
if((fwrite((void *)&rspi_pre[image_i], sizeof(struct resamp_info), 1, resampinfofile)) != 1 )
die("couldn't write to file", "resampinfo.bin");
fclose(resampinfofile);
////////////////////////////////////////////////////////////////////////////////////////////////////////////////
if (orb.points != NULL)
{
@ -198,12 +571,19 @@ void get_files(struct PRM *prm, FILE **rawfile, FILE **prmfile, char *prmfilenam
/* name and open output file for raw data (but input for later processing) */
/* if more than 1 set of output files, append an integer (beginning with 2) */
if (n == 0) {
sprintf(prm->input_file,"%s.raw", name);
//if (n == 0) {
// sprintf(prm->input_file,"%s.raw", name);
// sprintf(prmfilename,"%s.PRM", name);
//} else {
// sprintf(prm->input_file,"%s.raw.%d",name,n+1);
// sprintf(prmfilename,"%s.PRM.%d", name, n+1);
//}
if (n==0) {
sprintf(prmfilename,"%s.PRM", name);
sprintf(prm->input_file,"%s",name);
} else {
sprintf(prm->input_file,"%s.raw.%d",name,n+1);
sprintf(prmfilename,"%s.PRM.%d", name, n+1);
sprintf(prm->input_file,"%s.%d",name,n+1);
}
/* now open the files */
@ -212,4 +592,4 @@ void get_files(struct PRM *prm, FILE **rawfile, FILE **prmfile, char *prmfilenam
if ((*prmfile = fopen(prmfilename, "w")) == NULL) die ("couldn't open output PRM file \n",prmfilename);
}
/*------------------------------------------------------*/

View File

@ -6,11 +6,11 @@ Import('envSensorSrc1')
package = envSensorSrc1['PACKAGE']
project = envSensorSrc1['PROJECT']
install = envSensorSrc1['PRJ_LIB_DIR']
headerFiles = ['data_ALOS.h','data_ALOSE.h','image_sio.h','orbit_ALOS.h','sarleader_ALOS.h','sarleader_fdr.h','siocomplex.h']
sourceFiles = ['ALOSE_orbits_utils.c','ALOS_ldr_orbit.c','ALOS_pre_process.c','calc_dop.c','hermite_c.c','init_from_PRM.c',
'interpolate_ALOS_orbit.c','null_sio_struct.c','parse_ALOS_commands.c','polyfit.c','read_ALOSE_data.c',
headerFiles = ['data_ALOS.h','data_ALOSE.h','image_sio.h','orbit_ALOS.h','sarleader_ALOS.h','sarleader_fdr.h','siocomplex.h', 'resamp.h']
sourceFiles = ['ALOSE_orbits_utils.c','ALOS_ldr_orbit.c','ALOS_pre_process.c','calc_dop.c','get_sio_struct.c','hermite_c.c','init_from_PRM.c',
'interpolate_ALOS_orbit.c','null_sio_struct.c','parse_ALOS_commands.c','polyfit.c','put_sio_struct.c','read_ALOSE_data.c',
'read_ALOS_data.c','read_ALOS_sarleader.c','roi_utils.c','set_ALOS_defaults.c','siocomplex.c',
'swap_ALOS_data_info.c','utils.c','write_ALOS_prm.c','readOrbitPulse.f','readOrbitPulseState.f','readOrbitPulseSetState.f']
'swap_ALOS_data_info.c','utils.c','write_ALOS_prm.c','readOrbitPulse.f','readOrbitPulseState.f','readOrbitPulseSetState.f', 'lib_array.c', 'lib_cpx.c', 'lib_file.c', 'lib_func.c', 'resamp_azimuth.c']
lib = envSensorSrc1.Library(target = 'alos', source = sourceFiles)
envSensorSrc1.Install(install,lib)
envSensorSrc1.Alias('install',install)

View File

@ -12,6 +12,9 @@
* Date: *
* *****************************************************************************/
//program updated to handle PRF change issue of ALOS-1
//Cunren Liang, December 2019
#include "image_sio.h"
#include "lib_functions.h"
#include "siocomplex.h"
@ -23,8 +26,8 @@ void calc_dop(struct PRM *prm)
long n;
float *xr, *ac, *sg;
double sumd;
fcomplex *ai, *bi, *ab;
fcomplex ctmp;
fcomplex_sio *ai, *bi, *ab;
fcomplex_sio ctmp;
FILE *fin;
fprintf(stderr,".... calculating doppler for %s\n",prm->input_file);
@ -40,9 +43,15 @@ void calc_dop(struct PRM *prm)
ac = (float *) malloc(n*sizeof(float));
sg = (float *) malloc(n*sizeof(float));
ai = (fcomplex *) malloc(n*sizeof(fcomplex));
bi = (fcomplex *) malloc(n*sizeof(fcomplex));
ab = (fcomplex *) malloc(2*n*sizeof(fcomplex));
ai = (fcomplex_sio *) malloc(n*sizeof(fcomplex_sio));
bi = (fcomplex_sio *) malloc(n*sizeof(fcomplex_sio));
ab = (fcomplex_sio *) malloc(2*n*sizeof(fcomplex_sio));
for(i = 0; i< n;i++){
ab[i].r = 0;
ab[i].i = 0;
}
/* read a line of data from fin (input file, chars) to ai (complex floats) */
fread(indata, sizeof(unsigned char), prm->bytes_per_line, fin);
@ -52,7 +61,7 @@ void calc_dop(struct PRM *prm)
/* inefficient; could put loops inside each other */
for (i=prm->first_line; i<prm->num_lines-1; i++){
if (i/2000 == i/2000.0) fprintf(stderr," Working on line %d \n",i);
//if (i/2000 == i/2000.0) fprintf(stderr," Working on line %d \n",i);
fread(indata, sizeof(unsigned char), prm->bytes_per_line, fin);
@ -87,9 +96,10 @@ void calc_dop(struct PRM *prm)
free(xr); free(ac); free(sg);
free(ai); free(bi); free(ab);
free(indata);
fprintf(stderr,"done\n");
}
/*---------------------------------------------------*/
void read_data(fcomplex *data, unsigned char *indata, int i, struct PRM *prm)
void read_data(fcomplex_sio *data, unsigned char *indata, int i, struct PRM *prm)
{
int ii ;

View File

@ -0,0 +1,201 @@
/*--------------------------------------------------------------------*/
/*
Read parameters into PRM structure from PRM file
Based on get_params by Evelyn J. Price
Modified by RJM
*/
/*--------------------------------------------------------------------*/
//program updated to handle PRF change issue of ALOS-1
//Cunren Liang, December 2019
#include "image_sio.h"
#include "lib_functions.h"
/*
void get_sio_struct(FILE *, struct PRM *);
void get_string(char *, char *, char *, char *);
void get_int(char *, char *, char *, int *);
void get_double(char *, char *, char *, double *);
*/
void get_sio_struct(FILE *fh, struct PRM *s) {
char name[256], value[256];
debug = 0;
if (debug) {
fprintf(stderr, "get_sio_struct:\n");
fprintf(stderr, "PRMname (PRM value) interpreted value\n");
}
while (fscanf(fh, "%s = %s \n", name, value) != EOF) {
/* strings */
if (strcmp(name, "input_file") == 0)
get_string(name, "input_file", value, s->input_file);
if (strcmp(name, "led_file") == 0)
get_string(name, "led_file", value, s->led_file);
if (strcmp(name, "out_amp_file") == 0)
get_string(name, "out_amp_file", value, s->out_amp_file);
if (strcmp(name, "out_data_file") == 0)
get_string(name, "out_data_file", value, s->out_data_file);
if (strcmp(name, "scnd_rng_mig") == 0)
get_string(name, "scnd_rng_mig", value, s->srm);
if (strcmp(name, "deskew") == 0)
get_string(name, "deskew", value, s->deskew);
if (strcmp(name, "Flip_iq") == 0)
get_string(name, "Flip_iq", value, s->iqflip);
if (strcmp(name, "offset_video") == 0)
get_string(name, "offset_video", value, s->offset_video);
if (strcmp(name, "ref_file") == 0)
get_string(name, "ref_file", value, s->ref_file);
if (strcmp(name, "SLC_file") == 0)
get_string(name, "SLC_file", value, s->SLC_file);
if (strcmp(name, "orbdir") == 0)
get_string(name, "orbdir", value, s->orbdir);
//if (strcmp(name, "lookdir") == 0)
// get_string(name, "lookdir", value, s->lookdir);
if (strcmp(name, "date") == 0)
get_string(name, "date", value, s->date);
/* integers */
if (strcmp(name, "nrows") == 0)
get_int(name, "nrows", value, &s->nrows);
if (strcmp(name, "num_lines") == 0)
get_int(name, "num_lines", value, &s->num_lines);
if (strcmp(name, "bytes_per_line") == 0)
get_int(name, "bytes_per_line", value, &s->bytes_per_line);
if (strcmp(name, "good_bytes_per_line") == 0)
get_int(name, "good_bytes_per_line", value, &s->good_bytes);
if (strcmp(name, "first_line") == 0)
get_int(name, "first_line", value, &s->first_line);
if (strcmp(name, "num_patches") == 0)
get_int(name, "num_patches", value, &s->num_patches);
if (strcmp(name, "first_sample") == 0)
get_int(name, "first_sample", value, &s->first_sample);
if (strcmp(name, "num_valid_az") == 0)
get_int(name, "num_valid_az", value, &s->num_valid_az);
if (strcmp(name, "SC_identity") == 0)
get_int(name, "SC_identity", value, &s->SC_identity);
if (strcmp(name, "chirp_ext") == 0)
get_int(name, "chirp_ext", value, &s->chirp_ext);
if (strcmp(name, "st_rng_bin") == 0)
get_int(name, "st_rng_bin", value, &s->st_rng_bin);
if (strcmp(name, "num_rng_bins") == 0)
get_int(name, "num_rng_bins", value, &s->num_rng_bins);
if (strcmp(name, "ref_identity") == 0)
get_int(name, "ref_identity", value, &s->ref_identity);
if (strcmp(name, "nlooks") == 0)
get_int(name, "nlooks", value, &s->nlooks);
if (strcmp(name, "rshift") == 0)
get_int(name, "rshift", value, &s->rshift);
if (strcmp(name, "ashift") == 0)
get_int(name, "ashift", value, &s->ashift);
/* backwards compatibility for xshift/rshift yshift/ashift */
if (strcmp(name, "xshift") == 0)
get_int(name, "rshift", value, &s->rshift);
if (strcmp(name, "yshift") == 0)
get_int(name, "ashift", value, &s->ashift);
if (strcmp(name, "SLC_format") == 0)
get_int(name, "SLC_format", value, &s->SLC_format);
/* doubles */
if (strcmp(name, "SC_clock_start") == 0)
get_double(name, "SC_clock_start", value, &s->SC_clock_start);
if (strcmp(name, "SC_clock_stop") == 0)
get_double(name, "SC_clock_stop", value, &s->SC_clock_stop);
if (strcmp(name, "icu_start") == 0)
get_double(name, "icu_start", value, &s->icu_start);
//if (strcmp(name, "clock_start") == 0)
// get_double(name, "clock_start", value, &s->clock_start);
//if (strcmp(name, "clock_stop") == 0)
// get_double(name, "clock_stop", value, &s->clock_stop);
if (strcmp(name, "caltone") == 0)
get_double(name, "caltone", value, &s->caltone);
if (strcmp(name, "earth_radius") == 0)
get_double(name, "earth_radius", value, &s->RE);
if (strcmp(name, "equatorial_radius") == 0)
get_double(name, "equatorial_radius", value, &s->ra);
if (strcmp(name, "polar_radius") == 0)
get_double(name, "polar_radius", value, &s->rc);
if (strcmp(name, "SC_vel") == 0)
get_double(name, "SC_vel", value, &s->vel);
if (strcmp(name, "SC_height") == 0)
get_double(name, "SC_height", value, &s->ht);
if (strcmp(name, "SC_height_start") == 0)
get_double(name, "SC_height_start", value, &s->ht_start);
if (strcmp(name, "SC_height_end") == 0)
get_double(name, "SC_height_end", value, &s->ht_end);
if (strcmp(name, "near_range") == 0)
get_double(name, "near_range", value, &s->near_range);
if (strcmp(name, "PRF") == 0)
get_double(name, "PRF", value, &s->prf);
if (strcmp(name, "I_mean") == 0)
get_double(name, "I_mean", value, &s->xmi);
if (strcmp(name, "Q_mean") == 0)
get_double(name, "Q_mean", value, &s->xmq);
if (strcmp(name, "az_res") == 0)
get_double(name, "az_res", value, &s->az_res);
if (strcmp(name, "rng_samp_rate") == 0)
get_double(name, "rng_samp_rate", value, &s->fs);
if (strcmp(name, "chirp_slope") == 0)
get_double(name, "chirp_slope", value, &s->chirp_slope);
if (strcmp(name, "pulse_dur") == 0)
get_double(name, "pulse_dur", value, &s->pulsedur);
if (strcmp(name, "radar_wavelength") == 0)
get_double(name, "radar_wavelength", value, &s->lambda);
if (strcmp(name, "rng_spec_wgt") == 0)
get_double(name, "rng_spec_wgt", value, &s->rhww);
if (strcmp(name, "rm_rng_band") == 0)
get_double(name, "rm_rng_band", value, &s->pctbw);
if (strcmp(name, "rm_az_band") == 0)
get_double(name, "rm_az_band", value, &s->pctbwaz);
if (strcmp(name, "fd1") == 0)
get_double(name, "fd1", value, &s->fd1);
if (strcmp(name, "fdd1") == 0)
get_double(name, "fdd1", value, &s->fdd1);
if (strcmp(name, "fddd1") == 0)
get_double(name, "fddd1", value, &s->fddd1);
if (strcmp(name, "sub_int_r") == 0)
get_double(name, "sub_int_r", value, &s->sub_int_r);
if (strcmp(name, "sub_int_a") == 0)
get_double(name, "sub_int_a", value, &s->sub_int_a);
if (strcmp(name, "stretch_r") == 0)
get_double(name, "stretch_r", value, &s->stretch_r);
if (strcmp(name, "stretch_a") == 0)
get_double(name, "stretch_a", value, &s->stretch_a);
if (strcmp(name, "a_stretch_r") == 0)
get_double(name, "a_stretch_r", value, &s->a_stretch_r);
if (strcmp(name, "a_stretch_a") == 0)
get_double(name, "a_stretch_a", value, &s->a_stretch_a);
if (strcmp(name, "baseline_start") == 0)
get_double(name, "baseline_start", value, &s->baseline_start);
if (strcmp(name, "alpha_start") == 0)
get_double(name, "alpha_start", value, &s->alpha_start);
if (strcmp(name, "baseline_end") == 0)
get_double(name, "baseline_end", value, &s->baseline_end);
if (strcmp(name, "alpha_end") == 0)
get_double(name, "alpha_end", value, &s->alpha_end);
//if (strcmp(name, "SLC_scale") == 0)
// get_double(name, "SLC_scale", value, &s->SLC_scale);
}
}
/*--------------------------------------------------------------------------------*/
void get_string(char *s1, char *name, char *value, char *s2) {
strcpy(s2, value);
if (debug == 1)
fprintf(stderr, " %s (%s) = %s\n", s1, name, value);
}
/*--------------------------------------------------------------------------------*/
void get_int(char *s1, char *name, char *value, int *iparam) {
*iparam = atoi(value);
if (debug == 1)
fprintf(stderr, " %s (%s) = %s (%d)\n", s1, name, value, *iparam);
}
/*--------------------------------------------------------------------------------*/
void get_double(char *s1, char *name, char *value, double *param) {
*param = atof(value);
if (debug == 1)
fprintf(stderr, " %s (%s) = %s (%lf)\n", s1, name, value, *param);
}
/*--------------------------------------------------------------------------------*/

View File

@ -1,3 +1,6 @@
//program updated to handle PRF change issue of ALOS-1
//Cunren Liang, December 2019
/* taken from soi.h */
#include <stdio.h>
#include <stdlib.h>
@ -33,9 +36,9 @@
#define NULL_DOUBLE -99999.9999
#define NULL_CHAR "XXXXXXXX"
typedef struct SCOMPLEX {short r,i;} scomplex;
typedef struct FCOMPLEX {float r,i;} fcomplex;
typedef struct DCOMPLEX {double r,i;} dcomplex;
typedef struct SCOMPLEX_SIO {short r,i;} scomplex_sio;
typedef struct FCOMPLEX_SIO {float r,i;} fcomplex_sio;
typedef struct DCOMPLEX_SIO {double r,i;} dcomplex_sio;
struct PRM {
char input_file[256];
@ -121,6 +124,20 @@ struct PRM {
double bpara; /* parallel baseline - added by RJM */
double bperp; /* perpendicular baseline - added by RJM */
};
struct resamp_info {
//we assume there are no more than 20 prfs per image
int nPRF; //number of prfs, start with 1
int frame_counter_start[20];
int frame_counter_end[20];
int num_lines[20];
int num_bins[20];
double prf[20];
double SC_clock_start[20]; /* YYDDD.DDDD */
double fd1[20];
double fdd1[20];
double fddd1[20];
char input_file[20][256]; //we assume there are no more than 256 characters in the file name
};
/*
offset_video off_vid
chirp_ext nextend

View File

@ -0,0 +1,575 @@
//////////////////////////////////////
// Cunren Liang, NASA JPL/Caltech
// Copyright 2017
//////////////////////////////////////
#include "resamp.h"
/****************************************************************/
/* allocating arrays */
/****************************************************************/
signed char *vector_char(long nl, long nh)
/* allocate a signed char vector with subscript range v[nl..nh] */
{
signed char *v;
v=(signed char *)malloc((size_t) ((nh-nl+1+NR_END)*sizeof(signed char)));
if (!v){
fprintf(stderr,"Error: cannot allocate 1-D vector\n");
exit(1);
}
return v-nl+NR_END;
}
void free_vector_char(signed char *v, long nl, long nh)
/* free a signed char vector allocated with vector() */
{
free((FREE_ARG) (v+nl-NR_END));
}
unsigned char *vector_unchar(long nl, long nh)
/* allocate a unsigned char vector with subscript range v[nl..nh] */
{
unsigned char *v;
v=(unsigned char *)malloc((size_t) ((nh-nl+1+NR_END)*sizeof(unsigned char)));
if (!v){
fprintf(stderr,"Error: cannot allocate 1-D vector\n");
exit(1);
}
return v-nl+NR_END;
}
void free_vector_unchar(unsigned char *v, long nl, long nh)
/* free a unsigned char vector allocated with vector() */
{
free((FREE_ARG) (v+nl-NR_END));
}
int *vector_int(long nl, long nh)
/* allocate an int vector with subscript range v[nl..nh] */
{
int *v;
v=(int *)malloc((size_t) ((nh-nl+1+NR_END)*sizeof(int)));
if (!v) nrerror("Error: cannot allocate vector_int()");
return v-nl+NR_END;
}
void free_vector_int(int *v, long nl, long nh)
/* free an int vector allocated with ivector() */
{
free((FREE_ARG) (v+nl-NR_END));
}
float *vector_float(long nl, long nh)
/* allocate a float vector with subscript range v[nl..nh] */
{
float *v;
v=(float *)malloc((size_t) ((nh-nl+1+NR_END)*sizeof(float)));
if (!v){
fprintf(stderr,"Error: cannot allocate 1-D vector\n");
exit(1);
}
return v-nl+NR_END;
}
void free_vector_float(float *v, long nl, long nh)
/* free a float vector allocated with vector() */
{
free((FREE_ARG) (v+nl-NR_END));
}
double *vector_double(long nl, long nh)
/* allocate a double vector with subscript range v[nl..nh] */
{
double *v;
v=(double *)malloc((size_t) ((nh-nl+1+NR_END)*sizeof(double)));
if (!v){
fprintf(stderr,"Error: cannot allocate 1-D vector\n");
exit(1);
}
return v-nl+NR_END;
}
void free_vector_double(double *v, long nl, long nh)
/* free a double vector allocated with vector() */
{
free((FREE_ARG) (v+nl-NR_END));
}
fcomplex *vector_fcomplex(long nl, long nh)
/* allocate a fcomplex vector with subscript range v[nl..nh] */
{
fcomplex *v;
v=(fcomplex *)malloc((size_t) ((nh-nl+1+NR_END)*sizeof(fcomplex)));
if (!v) nrerror("cannot allocate fcvector()");
return v-nl+NR_END;
}
void free_vector_fcomplex(fcomplex *v, long nl, long nh)
/* free a fcomplex vector allocated with fcvector() */
{
free((FREE_ARG) (v+nl-NR_END));
}
signed char **matrix_char(long nrl, long nrh, long ncl, long nch)
/* allocate a signed char matrix with subscript range m[nrl..nrh][ncl..nch] */
{
long i, nrow=nrh-nrl+1,ncol=nch-ncl+1;
signed char **m;
/* allocate pointers to rows */
m=(signed char **) malloc((size_t)((nrow+NR_END)*sizeof(signed char*)));
if (!m) nrerror("Error: cannot allocate vector2d_float()");
m += NR_END;
m -= nrl;
/* allocate rows and set pointers to them */
m[nrl]=(signed char *) malloc((size_t)((nrow*ncol+NR_END)*sizeof(signed char)));
if (!m[nrl]) nrerror("Error: cannot allocate vector2d_float()");
m[nrl] += NR_END;
m[nrl] -= ncl;
for(i=nrl+1;i<=nrh;i++) m[i]=m[i-1]+ncol;
/* return pointer to array of pointers to rows */
return m;
}
void free_matrix_char(signed char **m, long nrl, long nrh, long ncl, long nch)
/* free a signed char matrix allocated by matrix() */
{
free((FREE_ARG) (m[nrl]+ncl-NR_END));
free((FREE_ARG) (m+nrl-NR_END));
}
unsigned char **matrix_unchar(long nrl, long nrh, long ncl, long nch)
/* allocate a unsigned char matrix with subscript range m[nrl..nrh][ncl..nch] */
{
long i, nrow=nrh-nrl+1,ncol=nch-ncl+1;
unsigned char **m;
/* allocate pointers to rows */
m=(unsigned char **) malloc((size_t)((nrow+NR_END)*sizeof(unsigned char*)));
if (!m) nrerror("Error: cannot allocate vector2d_float()");
m += NR_END;
m -= nrl;
/* allocate rows and set pointers to them */
m[nrl]=(unsigned char *) malloc((size_t)((nrow*ncol+NR_END)*sizeof(unsigned char)));
if (!m[nrl]) nrerror("Error: cannot allocate vector2d_float()");
m[nrl] += NR_END;
m[nrl] -= ncl;
for(i=nrl+1;i<=nrh;i++) m[i]=m[i-1]+ncol;
/* return pointer to array of pointers to rows */
return m;
}
void free_matrix_unchar(unsigned char **m, long nrl, long nrh, long ncl, long nch)
/* free a unsigned char matrix allocated by matrix() */
{
free((FREE_ARG) (m[nrl]+ncl-NR_END));
free((FREE_ARG) (m+nrl-NR_END));
}
float **matrix_float(long nrl, long nrh, long ncl, long nch)
/* allocate a float matrix with subscript range m[nrl..nrh][ncl..nch] */
{
long i, nrow=nrh-nrl+1,ncol=nch-ncl+1;
float **m;
/* allocate pointers to rows */
m=(float **) malloc((size_t)((nrow+NR_END)*sizeof(float*)));
if (!m) nrerror("Error: cannot allocate vector2d_float()");
m += NR_END;
m -= nrl;
/* allocate rows and set pointers to them */
m[nrl]=(float *) malloc((size_t)((nrow*ncol+NR_END)*sizeof(float)));
if (!m[nrl]) nrerror("Error: cannot allocate vector2d_float()");
m[nrl] += NR_END;
m[nrl] -= ncl;
for(i=nrl+1;i<=nrh;i++) m[i]=m[i-1]+ncol;
/* return pointer to array of pointers to rows */
return m;
}
void free_matrix_float(float **m, long nrl, long nrh, long ncl, long nch)
/* free a float matrix allocated by matrix() */
{
free((FREE_ARG) (m[nrl]+ncl-NR_END));
free((FREE_ARG) (m+nrl-NR_END));
}
double **matrix_double(long nrl, long nrh, long ncl, long nch)
/* allocate a double matrix with subscript range m[nrl..nrh][ncl..nch] */
{
long i, nrow=nrh-nrl+1,ncol=nch-ncl+1;
double **m;
/* allocate pointers to rows */
m=(double **) malloc((size_t)((nrow+NR_END)*sizeof(double*)));
if (!m) nrerror("Error: cannot allocate vector2d_double()");
m += NR_END;
m -= nrl;
/* allocate rows and set pointers to them */
m[nrl]=(double *) malloc((size_t)((nrow*ncol+NR_END)*sizeof(double)));
if (!m[nrl]) nrerror("Error: cannot allocate vector2d_double()");
m[nrl] += NR_END;
m[nrl] -= ncl;
for(i=nrl+1;i<=nrh;i++) m[i]=m[i-1]+ncol;
/* return pointer to array of pointers to rows */
return m;
}
void free_matrix_double(double **m, long nrl, long nrh, long ncl, long nch)
/* free a double matrix allocated by matrix() */
{
free((FREE_ARG) (m[nrl]+ncl-NR_END));
free((FREE_ARG) (m+nrl-NR_END));
}
/****************************************************************/
/* allocating C-style arrays */
/****************************************************************/
FILE **array1d_FILE(long nc){
FILE **fv;
fv = (FILE **)malloc(nc * sizeof(FILE *));
if(!fv){
fprintf(stderr,"Error: cannot allocate 1-D FILE array\n");
exit(1);
}
return fv;
}
void free_array1d_FILE(FILE **fv){
free(fv);
}
signed char *array1d_char(long nc){
signed char *fv;
fv = (signed char*) malloc(nc * sizeof(signed char));
if(!fv){
fprintf(stderr,"Error: cannot allocate 1-D signed char vector\n");
exit(1);
}
return fv;
}
void free_array1d_char(signed char *fv){
free(fv);
}
unsigned char *array1d_unchar(long nc){
unsigned char *fv;
fv = (unsigned char*) malloc(nc * sizeof(unsigned char));
if(!fv){
fprintf(stderr,"Error: cannot allocate 1-D unsigned char vector\n");
exit(1);
}
return fv;
}
void free_array1d_unchar(unsigned char *fv){
free(fv);
}
int *array1d_int(long nc){
int *fv;
fv = (int*) malloc(nc * sizeof(int));
if(!fv){
fprintf(stderr,"Error: cannot allocate 1-D int array\n");
exit(1);
}
return fv;
}
void free_array1d_int(int *fv){
free(fv);
}
float *array1d_float(long nc){
float *fv;
fv = (float*) malloc(nc * sizeof(float));
if(!fv){
fprintf(stderr,"Error: cannot allocate 1-D float vector\n");
exit(1);
}
return fv;
}
void free_array1d_float(float *fv){
free(fv);
}
double *array1d_double(long nc){
double *fv;
fv = (double*) malloc(nc * sizeof(double));
if(!fv){
fprintf(stderr,"Error: cannot allocate 1-D double vector\n");
exit(1);
}
return fv;
}
void free_array1d_double(double *fv){
free(fv);
}
fcomplex *array1d_fcomplex(long nc){
fcomplex *fcv;
fcv = (fcomplex*) malloc(nc * sizeof(fcomplex));
if(!fcv){
fprintf(stderr,"Error: cannot allocate 1-D float complex vector\n");
exit(1);
}
return fcv;
}
void free_array1d_fcomplex(fcomplex *fcv){
free(fcv);
}
dcomplex *array1d_dcomplex(long nc){
dcomplex *fcv;
fcv = (dcomplex*) malloc(nc * sizeof(dcomplex));
if(!fcv){
fprintf(stderr,"Error: cannot allocate 1-D double complex vector\n");
exit(1);
}
return fcv;
}
void free_array1d_dcomplex(dcomplex *fcv){
free(fcv);
}
signed char **array2d_char(long nl, long nc){
/* allocate a signed char 2-D matrix */
signed char **m;
int i;
/* allocate pointers to rows */
m = (signed char **) malloc(nl * sizeof(signed char *));
if(!m){
fprintf(stderr,"Error: cannot allocate 2-D matrix\n");
exit(1);
}
/* allocate rows */
m[0] = (signed char*) malloc(nl * nc * sizeof(signed char));
if(!m[0]){
fprintf(stderr,"Error: cannot allocate 2-D matrix\n");
exit(1);
}
/* set pointers */
for(i = 1; i < nl; i++){
m[i] = m[i-1] + nc;
}
return m;
}
void free_array2d_char(signed char **m){
/* free a signed char matrix allocated by farray2d() */
free(m[0]);
free(m);
}
unsigned char **array2d_unchar(long nl, long nc){
/* allocate a unsigned char 2-D matrix */
unsigned char **m;
int i;
/* allocate pointers to rows */
m = (unsigned char **) malloc(nl * sizeof(unsigned char *));
if(!m){
fprintf(stderr,"Error: cannot allocate 2-D matrix\n");
exit(1);
}
/* allocate rows */
m[0] = (unsigned char*) malloc(nl * nc * sizeof(unsigned char));
if(!m[0]){
fprintf(stderr,"Error: cannot allocate 2-D matrix\n");
exit(1);
}
/* set pointers */
for(i = 1; i < nl; i++){
m[i] = m[i-1] + nc;
}
return m;
}
void free_array2d_unchar(unsigned char **m){
/* free a signed unchar matrix allocated by farray2d() */
free(m[0]);
free(m);
}
float **array2d_float(long nl, long nc){
/* allocate a float 2-D matrix */
float **m;
int i;
/* allocate pointers to rows */
m = (float **) malloc(nl * sizeof(float *));
if(!m){
fprintf(stderr,"Error: cannot allocate 2-D matrix\n");
exit(1);
}
/* allocate rows */
m[0] = (float*) malloc(nl * nc * sizeof(float));
if(!m[0]){
fprintf(stderr,"Error: cannot allocate 2-D matrix\n");
exit(1);
}
/* set pointers */
for(i = 1; i < nl; i++){
m[i] = m[i-1] + nc;
}
return m;
}
void free_array2d_float(float **m){
/* free a float matrix allocated by farray2d() */
free(m[0]);
free(m);
}
double **array2d_double(long nl, long nc){
/* allocate a double 2-D matrix */
double **m;
int i;
/* allocate pointers to rows */
m = (double **) malloc(nl * sizeof(double *));
if(!m){
fprintf(stderr,"Error: cannot allocate 2-D matrix\n");
exit(1);
}
/* allocate rows */
m[0] = (double*) malloc(nl * nc * sizeof(double));
if(!m[0]){
fprintf(stderr,"Error: cannot allocate 2-D matrix\n");
exit(1);
}
/* set pointers */
for(i = 1; i < nl; i++){
m[i] = m[i-1] + nc;
}
return m;
}
void free_array2d_double(double **m){
/* free a double matrix allocated by farray2d() */
free(m[0]);
free(m);
}
fcomplex **array2d_fcomplex(long nl, long nc){
/* allocate a fcomplex 2-D matrix */
fcomplex **m;
int i;
/* allocate pointers to rows */
m = (fcomplex **) malloc(nl * sizeof(fcomplex *));
if(!m){
fprintf(stderr,"Error: cannot allocate 2-D matrix\n");
exit(1);
}
/* allocate rows */
m[0] = (fcomplex*) malloc(nl * nc * sizeof(fcomplex));
if(!m[0]){
fprintf(stderr,"Error: cannot allocate 2-D matrix\n");
exit(1);
}
/* set pointers */
for(i = 1; i < nl; i++){
m[i] = m[i-1] + nc;
}
return m;
}
void free_array2d_fcomplex(fcomplex **m){
/* free a fcomplex matrix allocated by fcarray2d() */
free(m[0]);
free(m);
}
/****************************************************************/
/* handling error */
/****************************************************************/
void nrerror(char error_text[])
/* Numerical Recipes standard error handler */
{
fprintf(stderr,"Numerical Recipes run-time error...\n");
fprintf(stderr,"%s\n",error_text);
fprintf(stderr,"...now exiting to system...\n");
exit(1);
}

View File

@ -0,0 +1,72 @@
//////////////////////////////////////
// Cunren Liang, NASA JPL/Caltech
// Copyright 2017
//////////////////////////////////////
#include "resamp.h"
// complex operations
fcomplex cmul(fcomplex a, fcomplex b)
{
fcomplex c;
c.re=a.re*b.re-a.im*b.im;
c.im=a.im*b.re+a.re*b.im;
return c;
}
fcomplex cconj(fcomplex z)
{
fcomplex c;
c.re=z.re;
c.im = -z.im;
return c;
}
fcomplex cadd(fcomplex a, fcomplex b)
{
fcomplex c;
c.re=a.re+b.re;
c.im=a.im+b.im;
return c;
}
float xcabs(fcomplex z)
{
float x,y,ans,temp;
x=fabs(z.re);
y=fabs(z.im);
if (x == 0.0)
ans=y;
else if (y == 0.0)
ans=x;
else if (x > y) {
temp=y/x;
ans=x*sqrt(1.0+temp*temp);
} else {
temp=x/y;
ans=y*sqrt(1.0+temp*temp);
}
return ans;
}
float cphs(fcomplex z){
float ans;
if(z.re == 0.0 && z.im == 0.0)
ans = 0.0;
else
ans = atan2(z.im, z.re);
return ans;
//it seems that there is no need to add the if clause
//do a test:
// printf("%12.4f, %12.4f, %12.4f, %12.4f, %12.4f\n", \
// atan2(0.0, 1.0), atan2(1.0, 0.0), atan2(0.0, -1.0), atan2(-1.0, 0.0), atan2(0.0, 0.0));
//output:
// 0.0000, 1.5708, 3.1416, -1.5708, 0.0000
}

View File

@ -0,0 +1,43 @@
//////////////////////////////////////
// Cunren Liang, NASA JPL/Caltech
// Copyright 2017
//////////////////////////////////////
#include "resamp.h"
FILE *openfile(char *filename, char *pattern){
FILE *fp;
fp=fopen(filename, pattern);
if (fp==NULL){
fprintf(stderr,"Error: cannot open file: %s\n", filename);
exit(1);
}
return fp;
}
void readdata(void *data, size_t blocksize, FILE *fp){
if(fread(data, blocksize, 1, fp) != 1){
fprintf(stderr,"Error: cannot read data\n");
exit(1);
}
}
void writedata(void *data, size_t blocksize, FILE *fp){
if(fwrite(data, blocksize, 1, fp) != 1){
fprintf(stderr,"Error: cannot write data\n");
exit(1);
}
}
long file_length(FILE* fp, long width, long element_size){
long length;
fseeko(fp,0L,SEEK_END);
length = ftello(fp) / element_size / width;
rewind(fp);
return length;
}

View File

@ -0,0 +1,275 @@
//////////////////////////////////////
// Cunren Liang, NASA JPL/Caltech
// Copyright 2017
//////////////////////////////////////
#include "resamp.h"
long next_pow2(long a){
long i;
long x;
x = 2;
while(x < a){
x *= 2;
}
return x;
}
void circ_shift(fcomplex *in, int na, int nc){
int i;
int ncm;
ncm = nc%na;
if(ncm < 0){
for(i = 0; i < abs(ncm); i++)
left_shift(in, na);
}
else if(ncm > 0){
for(i = 0; i < ncm; i++)
right_shift(in, na);
}
else{ //ncm == 0, no need to shift
i = 0;
}
}
void left_shift(fcomplex *in, int na){
int i;
fcomplex x;
if(na < 1){
fprintf(stderr, "Error: array size < 1\n\n");
exit(1);
}
else if(na > 1){
x.re = in[0].re;
x.im = in[0].im;
for(i = 0; i <= na - 2; i++){
in[i].re = in[i+1].re;
in[i].im = in[i+1].im;
}
in[na-1].re = x.re;
in[na-1].im = x.im;
}
else{ //na==1, no need to shift
i = 0;
}
}
void right_shift(fcomplex *in, int na){
int i;
fcomplex x;
if(na < 1){
fprintf(stderr, "Error: array size < 1\n\n");
exit(1);
}
else if(na > 1){
x.re = in[na-1].re;
x.im = in[na-1].im;
for(i = na - 1; i >= 1; i--){
in[i].re = in[i-1].re;
in[i].im = in[i-1].im;
}
in[0].re = x.re;
in[0].im = x.im;
}
else{ //na==1, no need to shift
i = 0;
}
}
int roundfi(float a){
int b;
if(a > 0)
b = (int)(a + 0.5);
else if (a < 0)
b = (int)(a - 0.5);
else
b = a;
return b;
}
void sinc(int n, int m, float *coef){
int i;
int hmn;
hmn = n * m / 2;
for(i=-hmn; i<=hmn; i++){
if(i != 0){
coef[i] = sin(PI * i / m) / (PI * i / m);
//coef[i] = sin(pi * i / m) / (pi * i / m);
}
else{
coef[i] = 1.0;
}
}
}
//kaiser() is equivalent to kaiser2()
//it is created to just keep the same style of sinc().
void kaiser(int n, int m, float *coef, float beta){
int i;
int hmn;
float a;
hmn = n * m / 2;
for(i = -hmn; i <= hmn; i++){
a = 1.0 - 4.0 * i * i / (n * m) / (n * m);
coef[i] = bessi0(beta * sqrt(a)) / bessi0(beta);
}
}
void kaiser2(float beta, int n, float *coef){
int i;
int hn;
float a;
hn = (n - 1) / 2;
for(i = -hn; i<=hn; i++){
a = 1.0 - 4.0 * i * i / (n - 1.0) / (n - 1.0);
coef[i] = bessi0(beta * sqrt(a)) / bessi0(beta);
}
}
void bandpass_filter(float bw, float bc, int n, int nfft, int ncshift, float beta, fcomplex *filter){
int i;
float *kw;
int hn;
fcomplex bwx, bcx;
hn = (n-1)/2;
if(n > nfft){
fprintf(stderr, "Error: fft length too small!\n\n");
exit(1);
}
if(abs(ncshift) > nfft){
fprintf(stderr, "Error: fft length too small or shift too big!\n\n");
exit(1);
}
//set all the elements to zero
for(i = 0; i < nfft; i++){
filter[i].re = 0.0;
filter[i].im = 0.0;
}
//calculate kaiser window
kw = vector_float(-hn, hn);
kaiser2(beta, n, kw);
//calculate filter
for(i = -hn; i <= hn; i++){
bcx.re = cos(bc * 2.0 * PI * i);
bcx.im = sin(bc * 2.0 * PI * i);
if(i == 0){
bwx.re = 1.0;
bwx.im = 0.0;
}
else{
bwx.re = sin(bw * PI * i) / (bw * PI * i);
bwx.im = 0.0;
}
filter[i+hn] = cmul(bcx, bwx);
filter[i+hn].re = bw * kw[i] * filter[i+hn].re;
filter[i+hn].im = bw * kw[i] * filter[i+hn].im;
}
//circularly shift filter, we shift the filter to left.
ncshift = -abs(ncshift);
circ_shift(filter, nfft, ncshift);
free_vector_float(kw, -hn, hn);
}
float bessi0(float x)
{
float ax,ans;
double y;
if ((ax=fabs(x)) < 3.75) {
y=x/3.75;
y*=y;
ans=1.0+y*(3.5156229+y*(3.0899424+y*(1.2067492
+y*(0.2659732+y*(0.360768e-1+y*0.45813e-2)))));
} else {
y=3.75/ax;
ans=(exp(ax)/sqrt(ax))*(0.39894228+y*(0.1328592e-1
+y*(0.225319e-2+y*(-0.157565e-2+y*(0.916281e-2
+y*(-0.2057706e-1+y*(0.2635537e-1+y*(-0.1647633e-1
+y*0.392377e-2))))))));
}
return ans;
}
#define SWAP(a,b) tempr=(a);(a)=(b);(b)=tempr
void four1(float data[], unsigned long nn, int isign)
{
unsigned long n,mmax,m,j,istep,i;
double wtemp,wr,wpr,wpi,wi,theta;
float tempr,tempi;
n=nn << 1;
j=1;
for (i=1;i<n;i+=2) {
if (j > i) {
SWAP(data[j],data[i]);
SWAP(data[j+1],data[i+1]);
}
m=nn;
while (m >= 2 && j > m) {
j -= m;
m >>= 1;
}
j += m;
}
mmax=2;
while (n > mmax) {
istep=mmax << 1;
theta=isign*(6.28318530717959/mmax);
wtemp=sin(0.5*theta);
wpr = -2.0*wtemp*wtemp;
wpi=sin(theta);
wr=1.0;
wi=0.0;
for (m=1;m<mmax;m+=2) {
for (i=m;i<=n;i+=istep) {
j=i+mmax;
tempr=wr*data[j]-wi*data[j+1];
tempi=wr*data[j+1]+wi*data[j];
data[j]=data[i]-tempr;
data[j+1]=data[i+1]-tempi;
data[i] += tempr;
data[i+1] += tempi;
}
wr=(wtemp=wr)*wpr-wi*wpi+wr;
wi=wi*wpr+wtemp*wpi+wi;
}
mmax=istep;
}
}
#undef SWAP

View File

@ -1,3 +1,6 @@
//program updated to handle PRF change issue of ALOS-1
//Cunren Liang, December 2019
/* include files to define sarleader structure */
#include "data_ALOS.h"
#include "data_ALOSE.h"
@ -9,8 +12,8 @@
void ALOS_ldr_orbit(struct ALOS_ORB *, struct PRM *);
void calc_height_velocity(struct ALOS_ORB *, struct PRM *, double, double, double *, double *, double *, double *, double *);
void calc_dop(struct PRM *);
void cfft1d_(int *, fcomplex *, int *);
void read_data(fcomplex *, unsigned char *, int, struct PRM *);
void cfft1d_(int *, fcomplex_sio *, int *);
void read_data(fcomplex_sio *, unsigned char *, int, struct PRM *);
void null_sio_struct(struct PRM *);
void get_sio_struct(FILE *, struct PRM *);
void put_sio_struct(struct PRM, FILE *);

View File

@ -0,0 +1,171 @@
/********************************************************************************
* Creator: Rob Mellors and David T. Sandwell * (San Diego State University,
*Scripps Institution of Oceanography) * Date : 10/03/2007 *
********************************************************************************/
/********************************************************************************
* Modification history: * Date: *
* *****************************************************************************/
//program updated to handle PRF change issue of ALOS-1
//Cunren Liang, December 2019
#include "image_sio.h"
#include "lib_functions.h"
/*
#define OUTFILE stdout
*/
/***************************************************************************/
void put_sio_struct(struct PRM prm, FILE *OUTFILE) {
/* set by set_ALOS_defaults */
if (prm.num_valid_az != NULL_INT)
fprintf(OUTFILE, "num_valid_az = %d \n", prm.num_valid_az);
if (prm.nrows != NULL_INT)
fprintf(OUTFILE, "nrows = %d \n", prm.nrows);
if (prm.first_line != NULL_INT)
fprintf(OUTFILE, "first_line = %d \n", prm.first_line);
if (strncmp(prm.deskew, NULL_CHAR, 8) != 0)
fprintf(OUTFILE, "deskew = %s \n", prm.deskew);
if (prm.caltone != NULL_DOUBLE)
fprintf(OUTFILE, "caltone = %lf \n", prm.caltone);
if (prm.st_rng_bin != NULL_INT)
fprintf(OUTFILE, "st_rng_bin = %d \n", prm.st_rng_bin);
if (strncmp(prm.iqflip, NULL_CHAR, 8) != 0)
fprintf(OUTFILE, "Flip_iq = %s \n", prm.iqflip);
if (strncmp(prm.offset_video, NULL_CHAR, 8) != 0)
fprintf(OUTFILE, "offset_video = %s \n", prm.offset_video);
if (prm.az_res != NULL_DOUBLE)
fprintf(OUTFILE, "az_res = %lf \n", prm.az_res);
if (prm.nlooks != NULL_INT)
fprintf(OUTFILE, "nlooks = %d \n", prm.nlooks);
if (prm.chirp_ext != NULL_INT)
fprintf(OUTFILE, "chirp_ext = %d \n", prm.chirp_ext);
if (strncmp(prm.srm, NULL_CHAR, 8) != 0)
fprintf(OUTFILE, "scnd_rng_mig = %s \n", prm.srm);
if (prm.rhww != NULL_DOUBLE)
fprintf(OUTFILE, "rng_spec_wgt = %lf \n", prm.rhww);
if (prm.pctbw != NULL_DOUBLE)
fprintf(OUTFILE, "rm_rng_band = %lf \n", prm.pctbw);
if (prm.pctbwaz != NULL_DOUBLE)
fprintf(OUTFILE, "rm_az_band = %lf \n", prm.pctbwaz);
if (prm.rshift != NULL_INT)
fprintf(OUTFILE, "rshift = %d \n", prm.rshift);
if (prm.ashift != NULL_INT)
fprintf(OUTFILE, "ashift = %d \n", prm.ashift);
if (prm.stretch_a != NULL_DOUBLE)
fprintf(OUTFILE, "stretch_r = %lf \n", prm.stretch_r);
if (prm.stretch_a != NULL_DOUBLE)
fprintf(OUTFILE, "stretch_a = %lf \n", prm.stretch_a);
if (prm.a_stretch_r != NULL_DOUBLE)
fprintf(OUTFILE, "a_stretch_r = %lf \n", prm.a_stretch_r);
if (prm.a_stretch_a != NULL_DOUBLE)
fprintf(OUTFILE, "a_stretch_a = %lf \n", prm.a_stretch_a);
if (prm.first_sample != NULL_INT)
fprintf(OUTFILE, "first_sample = %d \n", prm.first_sample);
if (prm.SC_identity != NULL_INT)
fprintf(OUTFILE, "SC_identity = %d \n", prm.SC_identity);
if (prm.fs != NULL_DOUBLE)
fprintf(OUTFILE, "rng_samp_rate = %lf \n", prm.fs);
/* from read_ALOS_data */
if (strncmp(prm.input_file, NULL_CHAR, 8) != 0)
fprintf(OUTFILE, "input_file = %s \n", prm.input_file);
if (prm.num_rng_bins != NULL_INT)
fprintf(OUTFILE, "num_rng_bins = %d \n", prm.num_rng_bins);
if (prm.bytes_per_line != NULL_INT)
fprintf(OUTFILE, "bytes_per_line = %d \n", prm.bytes_per_line);
if (prm.good_bytes != NULL_INT)
fprintf(OUTFILE, "good_bytes_per_line = %d \n", prm.good_bytes);
if (prm.prf != NULL_DOUBLE)
fprintf(OUTFILE, "PRF = %lf \n", prm.prf);
if (prm.pulsedur != NULL_DOUBLE)
fprintf(OUTFILE, "pulse_dur = %e \n", prm.pulsedur);
if (prm.near_range != NULL_DOUBLE)
fprintf(OUTFILE, "near_range = %lf \n", prm.near_range);
if (prm.num_lines != NULL_INT)
fprintf(OUTFILE, "num_lines = %d \n", prm.num_lines);
if (prm.num_patches != NULL_INT)
fprintf(OUTFILE, "num_patches = %d \n", prm.num_patches);
if (prm.SC_clock_start != NULL_DOUBLE)
fprintf(OUTFILE, "SC_clock_start = %16.10lf \n", prm.SC_clock_start);
if (prm.SC_clock_stop != NULL_DOUBLE)
fprintf(OUTFILE, "SC_clock_stop = %16.10lf \n", prm.SC_clock_stop);
//if (prm.clock_start != NULL_DOUBLE)
// fprintf(OUTFILE, "clock_start = %16.12lf \n", prm.clock_start);
//if (prm.clock_stop != NULL_DOUBLE)
// fprintf(OUTFILE, "clock_stop = %16.12lf \n", prm.clock_stop);
if (strncmp(prm.led_file, NULL_CHAR, 8) != 0)
fprintf(OUTFILE, "led_file = %s \n", prm.led_file);
/* from read_ALOS_ldrfile */
if (strncmp(prm.date, NULL_CHAR, 8) != 0)
fprintf(OUTFILE, "date = %.6s \n", prm.date);
if (strncmp(prm.orbdir, NULL_CHAR, 8) != 0)
fprintf(OUTFILE, "orbdir = %.1s \n", prm.orbdir);
//if (strncmp(prm.lookdir, NULL_CHAR, 8) != 0)
// fprintf(OUTFILE, "lookdir = %.1s \n", prm.lookdir);
if (prm.lambda != NULL_DOUBLE)
fprintf(OUTFILE, "radar_wavelength = %lg \n", prm.lambda);
if (prm.chirp_slope != NULL_DOUBLE)
fprintf(OUTFILE, "chirp_slope = %lg \n", prm.chirp_slope);
if (prm.fs != NULL_DOUBLE)
fprintf(OUTFILE, "rng_samp_rate = %lf \n", prm.fs);
if (prm.xmi != NULL_DOUBLE)
fprintf(OUTFILE, "I_mean = %lg \n", prm.xmi);
if (prm.xmq != NULL_DOUBLE)
fprintf(OUTFILE, "Q_mean = %lg \n", prm.xmq);
if (prm.vel != NULL_DOUBLE)
fprintf(OUTFILE, "SC_vel = %lf \n", prm.vel);
if (prm.RE != NULL_DOUBLE)
fprintf(OUTFILE, "earth_radius = %lf \n", prm.RE);
if (prm.ra != NULL_DOUBLE)
fprintf(OUTFILE, "equatorial_radius = %lf \n", prm.ra);
if (prm.rc != NULL_DOUBLE)
fprintf(OUTFILE, "polar_radius = %lf \n", prm.rc);
if (prm.ht != NULL_DOUBLE)
fprintf(OUTFILE, "SC_height = %lf \n", prm.ht);
if (prm.ht_start != NULL_DOUBLE)
fprintf(OUTFILE, "SC_height_start = %lf \n", prm.ht_start);
if (prm.ht_end != NULL_DOUBLE)
fprintf(OUTFILE, "SC_height_end = %lf \n", prm.ht_end);
if (prm.fd1 != NULL_DOUBLE)
fprintf(OUTFILE, "fd1 = %lf \n", prm.fd1);
if (prm.fdd1 != NULL_DOUBLE)
fprintf(OUTFILE, "fdd1 = %12.8lf \n", prm.fdd1);
if (prm.fddd1 != NULL_DOUBLE)
fprintf(OUTFILE, "fddd1 = %lf \n", prm.fddd1);
/* from calc_baseline */
/*
if (prm.rshift != NULL_INT) fprintf(OUTFILE, "rshift = %d
\n",prm.rshift); if (prm.ashift != NULL_INT) fprintf(OUTFILE, "ashift =
%d\n",prm.ashift);
*/
if (prm.sub_int_r != NULL_DOUBLE)
fprintf(OUTFILE, "sub_int_r = %f \n", prm.sub_int_r);
if (prm.sub_int_a != NULL_DOUBLE)
fprintf(OUTFILE, "sub_int_a = %f \n", prm.sub_int_a);
if (prm.bpara != NULL_DOUBLE)
fprintf(OUTFILE, "B_parallel = %f \n", prm.bpara);
if (prm.bperp != NULL_DOUBLE)
fprintf(OUTFILE, "B_perpendicular = %f \n", prm.bperp);
if (prm.baseline_start != NULL_DOUBLE)
fprintf(OUTFILE, "baseline_start = %f \n", prm.baseline_start);
if (prm.alpha_start != NULL_DOUBLE)
fprintf(OUTFILE, "alpha_start = %f \n", prm.alpha_start);
if (prm.baseline_end != NULL_DOUBLE)
fprintf(OUTFILE, "baseline_end = %f \n", prm.baseline_end);
if (prm.alpha_end != NULL_DOUBLE)
fprintf(OUTFILE, "alpha_end = %f \n", prm.alpha_end);
/* from sarp */
if (strncmp(prm.SLC_file, NULL_CHAR, 8) != 0)
fprintf(OUTFILE, "SLC_file = %s \n", prm.SLC_file);
//if (strncmp(prm.dtype, NULL_CHAR, 8) != 0)
// fprintf(OUTFILE, "dtype = %s \n", prm.dtype);
//if (prm.SLC_scale != NULL_DOUBLE)
// fprintf(OUTFILE, "SLC_scale = %f \n", prm.SLC_scale);
}
/***************************************************************************/

View File

@ -22,6 +22,14 @@
* 15-Apr-2010 Replaced ALOS identifier with ALOSE Jeff Bytof *
**************************************************************************/
/********************************************************************************
This program has been upgraded to handle the ALOS-1 PRF change issue.
BUT HAS NOT BEEN TESTED YET!!!
Cunren Liang, 12-December-2019
California Institute of Technology, Pasadena, CA
*********************************************************************************/
/*
the data header information is read into the structure dfd
the line prefix information is read into sdr
@ -36,7 +44,8 @@ SC_clock_start SC_clock_stop
#include "image_sio.h"
#include "lib_functions.h"
#define ZERO_VALUE (char)(63 + rand() % 2)
#define clip127(A) (((A) > 127) ? 127 : (((A) < 0) ? 0 : A))
/*
#define znew (int) (z=36969*(z&65535)+(z>>16))
typedef unsigned long UL;
@ -54,11 +63,12 @@ int assign_sardata_params_ALOSE(struct PRM *, int, int *, int *);
void swap_ALOS_data_info(struct sardata_info_ALOSE *sdr);
void settable(unsigned long);
void print_params(struct PRM *prm);
int check_shift(struct PRM *, int *, int *, int *, int);
int check_shift(struct PRM *, int *, int *, int *, int, int);
int set_file_position(FILE *, long *, int);
int reset_params(struct PRM *prm, long *, int *, int *);
int fill_shift_data(int, int, int, int, int, char *, char *, FILE *);
int handle_prf_change_ALOSE(struct PRM *, FILE *, long *, int);
void change_dynamic_range(char *data, long length);
struct sardata_record r1;
struct sardata_descriptor_ALOSE dfd;
@ -74,11 +84,13 @@ struct sardata_info_ALOSE
SARDATA__WCS_ALOSE
SARDATA_RVL_ALOSE(SP)
*/
long read_ALOSE_data (FILE *imagefile, FILE *outfile, struct PRM *prm, long *byte_offset) {
long read_ALOSE_data (FILE *imagefile, FILE *outfile, struct PRM *prm, long *byte_offset, struct resamp_info *rspi, int nPRF) {
char *data_fbd, *data, *shift_data;
int record_length0; /* length of record read at start of file */
int record_length1; /* length of record read in file */
int start_sdr_rec_len = 0; /* sdr record length for fisrt record */
int slant_range_old = 0; /* slant range of previous record */
int line_suffix_size; /* number of bytes after data */
int data_length; /* bytes of data */
int k, n, m, ishift, shift, shift0;
@ -91,6 +103,13 @@ long read_ALOSE_data (FILE *imagefile, FILE *outfile, struct PRM *prm, long *byt
if (verbose) fprintf(stderr,".... reading header \n");
//here we still get sdr from the first data line no matter whether prf changes.
//this sdr is used to initialize record_length0 in assign_sardata_params, which
//is used at line 152 to check if record_length changed.
//I think we should get sdr from first prf-change data line for the output of prf-change file.
//Cunren Liang. 02-DEC-2019
/* read header information */
read_sardata_info_ALOSE(imagefile, prm, &header_size, &line_prefix_size);
if (verbose) fprintf(stderr,".... reading header %d %d\n", header_size, line_prefix_size);
@ -115,7 +134,7 @@ long read_ALOSE_data (FILE *imagefile, FILE *outfile, struct PRM *prm, long *byt
shift0 = 0;
n = 1;
m = 0;
m = 2;//first line sequence_number
/* read the rest of the file */
while ( (fread((void *) &sdr,sizeof(struct sardata_info_ALOSE), 1, imagefile)) == 1 ) {
@ -124,9 +143,26 @@ long read_ALOSE_data (FILE *imagefile, FILE *outfile, struct PRM *prm, long *byt
/* checks for little endian/ big endian */
if (swap) swap_ALOS_data_info(&sdr);
if (n == 2)
//rspi->frame_counter_start[nPRF] = sdr.frame_counter;
//unfortunately restec format does not have this info, so we are not able to adjust time
rspi->frame_counter_start[nPRF] = 0;
/* if this is partway through the file due to prf change, reset sequence, PRF, and near_range */
if (n == 2)
start_sdr_rec_len = sdr.record_length;
if ((*byte_offset > 0) && (n == 2)) reset_params(prm, byte_offset, &n, &m);
if (sdr.record_length != start_sdr_rec_len) {
printf(" ***** warning sdr.record_length error %d \n", sdr.record_length);
sdr.record_length = start_sdr_rec_len;
sdr.PRF = prm->prf;
sdr.slant_range = slant_range_old;
}
if (sdr.sequence_number != n) printf(" missing line: n, seq# %d %d \n", n, sdr.sequence_number);
/* check for changes in record_length and PRF */
@ -136,11 +172,15 @@ long read_ALOSE_data (FILE *imagefile, FILE *outfile, struct PRM *prm, long *byt
/* if prf changes, close file and set byte_offset */
if ((sdr.PRF) != prm->prf) {
handle_prf_change_ALOSE(prm, imagefile, byte_offset, n);
n-=1;
break;
}
//rspi->frame_counter_end[nPRF] = sdr.frame_counter;
//unfortunately restec format does not have this info, so we are not able to adjust time
rspi->frame_counter_end[nPRF] = 0;
/* check shift to see if it varies from beginning or from command line value */
check_shift(prm, &shift, &ishift, &shift0, record_length1);
check_shift(prm, &shift, &ishift, &shift0, record_length1, 1);
if ((verbose) && (n/2000.0 == n/2000)) {
fprintf(stderr," Working on line %d prf %f record length %d slant_range %d \n"
@ -151,6 +191,7 @@ long read_ALOSE_data (FILE *imagefile, FILE *outfile, struct PRM *prm, long *byt
if ( fread ((char *) data, record_length1, (size_t) 1, imagefile) != 1 ) break;
data_length = record_length1;
slant_range_old = sdr.slant_range;
/* write line header to output data */
/* PSA - turning off headers
@ -165,6 +206,7 @@ long read_ALOSE_data (FILE *imagefile, FILE *outfile, struct PRM *prm, long *byt
}
/* write fbd data */
if (shift == 0) {
change_dynamic_range(data_fbd, data_length/2);
fwrite((char *) data_fbd, data_length/2, 1, outfile);
} else if (shift != 0) {
fill_shift_data(shift, ishift, data_length/2, line_suffix_size, record_length1, data_fbd, shift_data, outfile);
@ -173,18 +215,27 @@ long read_ALOSE_data (FILE *imagefile, FILE *outfile, struct PRM *prm, long *byt
else {
/* write fbs data */
if (shift == 0) {
change_dynamic_range(data, data_length);
fwrite((char *) data, data_length, 1, outfile);
} else if (shift != 0) {
fill_shift_data(shift, ishift, data_length, line_suffix_size, record_length1, data, shift_data, outfile);
}
}
}
//we are not writing out line prefix data, need to correct these parameters
//as they are used in doppler computation.
prm->first_sample = 0;
prm->bytes_per_line -= line_prefix_size;
prm->good_bytes -= line_prefix_size;
//this is the sdr of the first prf-change data line, should seek back to get last sdr to be used here.
/* calculate end time */
prm->SC_clock_stop = get_clock_ALOSE(sdr, tbias);
/* m is non-zero only in the event of a prf change */
prm->num_lines = n - m - 1;
//not correct if PRF changes, so I updated it here.
prm->num_lines = n - m + 1;
prm->num_patches = (int)((1.0*n)/(1.0*prm->num_valid_az));
if (prm->num_lines == 0) prm->num_lines = 1;
@ -194,6 +245,15 @@ long read_ALOSE_data (FILE *imagefile, FILE *outfile, struct PRM *prm, long *byt
prm->prf = 1.e3/pri;
prm->xmi = 63.5;
prm->xmq = 63.5;
rspi->prf[nPRF] = prm->prf;
rspi->SC_clock_start[nPRF] = prm->SC_clock_start;
rspi->num_lines[nPRF] = prm->num_lines;
rspi->num_bins[nPRF] = prm->bytes_per_line/(2*sizeof(char));
if (verbose) print_params(prm);
free(data);
@ -321,14 +381,14 @@ double get_clock();
/***************************************************************************/
int handle_prf_change_ALOSE(struct PRM *prm, FILE *imagefile, long *byte_offset, int n)
{
prm->num_lines = n;
//prm->num_lines = n;
fseek(imagefile, -1*sizeof(struct sardata_info_ALOSE), SEEK_CUR);
*byte_offset = ftell(imagefile);
printf(" *** PRF changed from %lf to %lf at line %d (byte %ld)\n", (0.001*prm->prf),(0.001*sdr.PRF), n, *byte_offset);
printf(" end: PRF changed from %lf to %lf at line %d \n", (0.001*prm->prf),(0.001*sdr.PRF), n);
printf(" *** PRF changed from %lf to %lf at line %d (byte %ld)\n", (0.001*prm->prf),(0.001*sdr.PRF), n-1, *byte_offset);
// printf(" end: PRF changed from %lf to %lf at line %d \n", (0.001*prm->prf),(0.001*sdr.PRF), n);
return(EXIT_SUCCESS);
}

View File

@ -33,6 +33,13 @@
* 07/17/08 reformatted; added functions RJM *
***************************************************************************/
/********************************************************************************
This program has been upgraded to handle the ALOS-1 PRF change issue.
Cunren Liang, 12-December-2019
California Institute of Technology, Pasadena, CA
*********************************************************************************/
/*
the data header information is read into the structure dfd
the line prefix information is read into sdr
@ -47,7 +54,8 @@ SC_clock_start SC_clock_stop
#include "image_sio.h"
#include "lib_functions.h"
#define ZERO_VALUE (char)(63 + rand() % 2)
#define clip127(A) (((A) > 127) ? 127 : (((A) < 0) ? 0 : A))
#define znew (int) (z=36969*(z&65535)+(z>>16))
typedef unsigned long UL;
static UL z=362436069, t[256];
@ -61,21 +69,24 @@ void swap_ALOS_data_info(struct sardata_info *sdr);
long read_sardata_info(FILE *, struct PRM *, int *, int *);
void print_params(struct PRM *prm);
int assign_sardata_params(struct PRM *, int, int *, int *);
int check_shift(struct PRM *, int *, int *, int *, int);
int check_shift(struct PRM *, int *, int *, int *, int, int);
int set_file_position(FILE *, long *, int);
int reset_params(struct PRM *prm, long *, int *, int *);
int fill_shift_data(int, int, int, int, int, char *, char *, FILE *);
int handle_prf_change(struct PRM *, FILE *, long *, int);
void change_dynamic_range(char *data, long length);
struct sardata_record r1;
struct sardata_descriptor dfd;
struct sardata_info sdr;
long read_ALOS_data (FILE *imagefile, FILE *outfile, struct PRM *prm, long *byte_offset) {
long read_ALOS_data (FILE *imagefile, FILE *outfile, struct PRM *prm, long *byte_offset, struct resamp_info *rspi, int nPRF) {
char *data, *shift_data;
int record_length0; /* length of record read at start of file */
int record_length1; /* length of record read in file */
int start_sdr_rec_len = 0; /* sdr record length for fisrt record */
int slant_range_old = 0; /* slant range of previous record */
int line_suffix_size; /* number of bytes after data */
int data_length; /* bytes of data */
int n, m, ishift, shift, shift0, npatch_max;
@ -87,6 +98,13 @@ long read_ALOS_data (FILE *imagefile, FILE *outfile, struct PRM *prm, long *byte
if (debug) fprintf(stderr,".... reading header \n");
//here we still get sdr from the first data line no matter whether prf changes.
//this sdr is used to initialize record_length0 in assign_sardata_params, which
//is used at line 152 to check if record_length changed.
//I think we should get sdr from first prf-change data line for the output of prf-change file.
//Cunren Liang. 02-DEC-2019
/* read header information */
read_sardata_info(imagefile, prm, &header_size, &line_prefix_size);
@ -111,7 +129,7 @@ long read_ALOS_data (FILE *imagefile, FILE *outfile, struct PRM *prm, long *byte
shift0 = 0;
n = 1;
m = 0;
m = 2;//first line sequence_number
/* read the rest of the file */
while ( (fread((void *) &sdr,sizeof(struct sardata_info), 1, imagefile)) == 1 ) {
@ -120,10 +138,29 @@ long read_ALOS_data (FILE *imagefile, FILE *outfile, struct PRM *prm, long *byte
/* checks for little endian/ big endian */
if (swap) swap_ALOS_data_info(&sdr);
/* if this is partway through the file due to prf change, reset sequence, PRF, and near_range */
if ((*byte_offset > 0) && (n == 2)) reset_params(prm, byte_offset, &n, &m);
if (sdr.sequence_number != n) printf(" missing line: n, seq# %d %d \n", n, sdr.sequence_number);
if (n == 2)
rspi->frame_counter_start[nPRF] = sdr.frame_counter;
/* if this is partway through the file due to prf change, reset sequence,
* PRF, and near_range */
if (n == 2)
start_sdr_rec_len = sdr.record_length;
if ((*byte_offset > 0) && (n == 2))
reset_params(prm, byte_offset, &n, &m);
if (sdr.record_length != start_sdr_rec_len) {
printf(" ***** warning sdr.record_length error %d \n", sdr.record_length);
sdr.record_length = start_sdr_rec_len;
sdr.PRF = prm->prf;
sdr.slant_range = slant_range_old;
}
if (sdr.sequence_number != n)
printf(" missing line: n, seq# %d %d \n", n, sdr.sequence_number);
/* check for changes in record_length and PRF */
record_length1 = sdr.record_length - line_prefix_size;
@ -132,11 +169,13 @@ long read_ALOS_data (FILE *imagefile, FILE *outfile, struct PRM *prm, long *byte
/* if prf changes, close file and set byte_offset */
if ((sdr.PRF) != prm->prf) {
handle_prf_change(prm, imagefile, byte_offset, n);
n-=1;
break;
}
rspi->frame_counter_end[nPRF] = sdr.frame_counter;
/* check shift to see if it varies from beginning or from command line value */
check_shift(prm, &shift, &ishift, &shift0, record_length1);
check_shift(prm, &shift, &ishift, &shift0, record_length1, 0);
if ((verbose) && (n/2000.0 == n/2000)) {
fprintf(stderr," Working on line %d prf %f record length %d slant_range %d \n"
@ -147,27 +186,37 @@ long read_ALOS_data (FILE *imagefile, FILE *outfile, struct PRM *prm, long *byte
if ( fread ((char *) data, record_length1, (size_t) 1, imagefile) != 1 ) break;
data_length = record_length1;
slant_range_old = sdr.slant_range;
/* write line header to output data */
/* PSA turning off headers
fwrite((void *) &sdr, line_prefix_size, 1, outfile); */
//header is not written to output
//fwrite((void *) &sdr, line_prefix_size, 1, outfile);
/* write data */
if (shift == 0) {
change_dynamic_range(data, data_length);
fwrite((char *) data, data_length, 1, outfile);
/* if data is shifted, fill in with data values of NULL_DATA at start or end*/
} else if (shift != 0) {
fill_shift_data(shift, ishift, data_length, line_suffix_size, record_length1, data, shift_data, outfile);
}
}
//we are not writing out line prefix data, need to correct these parameters
//as they are used in doppler computation.
prm->first_sample = 0;
prm->bytes_per_line -= line_prefix_size;
prm->good_bytes -= line_prefix_size;
/* calculate end time and fix prf */
prm->prf = 0.001*prm->prf;
//this is the sdr of the first prf-change data line, should seek back to get last sdr to be used here.
prm->SC_clock_stop = get_clock(sdr, tbias);
/* m is non-zero only in the event of a prf change */
prm->num_lines = n - m - 1;
//not correct if PRF changes, so I updated it here.
prm->num_lines = n - m + 1;
/* calculate the maximum number of patches and use that if the default is set to 1000 */
npatch_max = (int)((1.0*n)/(1.0*prm->num_valid_az));
@ -175,6 +224,15 @@ long read_ALOS_data (FILE *imagefile, FILE *outfile, struct PRM *prm, long *byte
if (prm->num_lines == 0) prm->num_lines = 1;
prm->xmi = 63.5;
prm->xmq = 63.5;
rspi->prf[nPRF] = prm->prf;
rspi->SC_clock_start[nPRF] = prm->SC_clock_start;
rspi->num_lines[nPRF] = prm->num_lines;
rspi->num_bins[nPRF] = prm->bytes_per_line/(2*sizeof(char));
if (verbose) print_params(prm);
free(data);
@ -293,7 +351,7 @@ double get_clock();
return(EXIT_SUCCESS);
}
/***************************************************************************/
int check_shift(struct PRM *prm, int *shift, int *ishift, int *shift0, int record_length1)
int check_shift(struct PRM *prm, int *shift, int *ishift, int *shift0, int record_length1, int ALOS_format)
{
*shift = 2*floor(0.5 + (sdr.slant_range - prm->near_range)/(0.5*SOL/prm->fs));
*ishift = abs(*shift);
@ -304,7 +362,13 @@ int check_shift(struct PRM *prm, int *shift, int *ishift, int *shift0, int recor
}
if(*shift != *shift0) {
printf(" near_range, shift = %d %d \n", sdr.slant_range, *shift);
if(ALOS_format==0)
printf(" near_range, shift = %d %d , at frame_counter: %d, line number: %d\n", sdr.slant_range, *shift, sdr.frame_counter, sdr.sequence_number-1);
if(ALOS_format==1)
printf(" near_range, shift = %d %d\n", sdr.slant_range, *shift);
*shift0 = *shift;
}
@ -324,21 +388,21 @@ int set_file_position(FILE *imagefile, long *byte_offset, int header_size)
return(EXIT_SUCCESS);
}
/***************************************************************************/
int reset_params(struct PRM *prm, long *byte_offset, int *n, int *m)
{
double get_clock();
int reset_params(struct PRM *prm, long *byte_offset, int *n, int *m) {
double get_clock();
prm->SC_clock_start = get_clock(sdr, tbias);
prm->prf = sdr.PRF;
prm->near_range = sdr.slant_range;
//comment out so that all data files with different prfs can be aligned at the same starting range
//prm->near_range = sdr.slant_range;
*n = sdr.sequence_number;
*m = *n;
*byte_offset = 0;
if (verbose) {
fprintf(stderr," new parameters: \n sequence number %d \n PRF %f\n near_range %lf\n",
*n, 0.001*prm->prf,prm->near_range);
}
return(EXIT_SUCCESS);
fprintf(stderr, " new parameters: \n sequence number %d \n PRF %f\n near_range %lf\n", *n, 0.001 * prm->prf,
prm->near_range);
}
return (EXIT_SUCCESS);
}
/***************************************************************************/
int fill_shift_data(int shift, int ishift, int data_length,
@ -359,6 +423,7 @@ int k;
}
/* write the shifted data out */
change_dynamic_range(shift_data, data_length);
fwrite((char *) shift_data, data_length, 1, outfile);
return(EXIT_SUCCESS);
@ -366,7 +431,7 @@ int k;
/***************************************************************************/
int handle_prf_change(struct PRM *prm, FILE *imagefile, long *byte_offset, int n)
{
prm->num_lines = n;
//prm->num_lines = n;
/* skip back to beginning of the line */
fseek(imagefile, -1*sizeof(struct sardata_info), SEEK_CUR);
@ -375,9 +440,34 @@ int handle_prf_change(struct PRM *prm, FILE *imagefile, long *byte_offset, int n
*byte_offset = ftell(imagefile);
/* tell the world */
printf(" *** PRF changed from %lf to %lf at line %d (byte %ld)\n", (0.001*prm->prf),(0.001*sdr.PRF), n, *byte_offset);
printf(" end: PRF changed from %lf to %lf at line %d \n", (0.001*prm->prf),(0.001*sdr.PRF), n);
printf(" *** PRF changed from %lf to %lf at line %d (byte %ld)\n", (0.001*prm->prf),(0.001*sdr.PRF), n-1, *byte_offset);
// printf(" end: PRF changed from %lf to %lf at line %d \n", (0.001*prm->prf),(0.001*sdr.PRF), n);
return(EXIT_SUCCESS);
}
/***************************************************************************/
void change_dynamic_range(char *data, long length){
long i;
for(i = 0; i < length; i++)
//THIS SHOULD NOT AFFECT DOPPLER COMPUTATION (SUCH AS IN calc_dop.c), BECAUSE
// 1. IQ BIAS IS REMOVED BEFORE COMPUTATION OF DOPPLER.
// 2. 2.0 WILL BE CANCELLED OUT IN atan2f().
// 3. actual computation results also verified this (even if there is a difference, it is about 0.* Hz)
//data[i] = (unsigned char)clip127(rintf(2. * (data[i] - 15.5) + 63.5));
data[i] = (unsigned char)clip127(rintf(2.0 * (data[i] - 15.5) + ZERO_VALUE));
}

View File

@ -0,0 +1,106 @@
//////////////////////////////////////
// Cunren Liang, NASA JPL/Caltech
// Copyright 2017
//////////////////////////////////////
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <math.h>
#define NR_END 1
#define FREE_ARG char*
#define PI 3.1415926535897932384626433832795028841971693993751058
typedef struct {
float re;
float im;
} fcomplex;
typedef struct {
double re;
double im;
} dcomplex;
//allocate arrays
signed char *vector_char(long nl, long nh);
void free_vector_char(signed char *v, long nl, long nh);
unsigned char *vector_unchar(long nl, long nh);
void free_vector_unchar(unsigned char *v, long nl, long nh);
int *vector_int(long nl, long nh);
void free_vector_int(int *v, long nl, long nh);
float *vector_float(long nl, long nh);
void free_vector_float(float *v, long nl, long nh);
double *vector_double(long nl, long nh);
void free_vector_double(double *v, long nl, long nh);
fcomplex *vector_fcomplex(long nl, long nh);
void free_vector_fcomplex(fcomplex *v, long nl, long nh);
signed char **matrix_char(long nrl, long nrh, long ncl, long nch);
void free_matrix_char(signed char **m, long nrl, long nrh, long ncl, long nch);
unsigned char **matrix_unchar(long nrl, long nrh, long ncl, long nch);
void free_matrix_unchar(unsigned char **m, long nrl, long nrh, long ncl, long nch);
float **matrix_float(long nrl, long nrh, long ncl, long nch);
void free_matrix_float(float **m, long nrl, long nrh, long ncl, long nch);
double **matrix_double(long nrl, long nrh, long ncl, long nch);
void free_matrix_double(double **m, long nrl, long nrh, long ncl, long nch);
//allocate C-style arrays
FILE **array1d_FILE(long nc);
void free_array1d_FILE(FILE **fv);
signed char *array1d_char(long nc);
void free_array1d_char(signed char *fv);
unsigned char *array1d_unchar(long nc);
void free_array1d_unchar(unsigned char *fv);
int *array1d_int(long nc);
void free_array1d_int(int *fv);
float *array1d_float(long nc);
void free_array1d_float(float *fv);
double *array1d_double(long nc);
void free_array1d_double(double *fv);
fcomplex *array1d_fcomplex(long nc);
void free_array1d_fcomplex(fcomplex *fcv);
dcomplex *array1d_dcomplex(long nc);
void free_array1d_dcomplex(dcomplex *fcv);
signed char **array2d_char(long nl, long nc);
void free_array2d_char(signed char **m);
unsigned char **array2d_unchar(long nl, long nc);
void free_array2d_unchar(unsigned char **m);
float **array2d_float(long nl, long nc);
void free_array2d_float(float **m);
double **array2d_double(long nl, long nc);
void free_array2d_double(double **m);
fcomplex **array2d_fcomplex(long nl, long nc);
void free_array2d_fcomplex(fcomplex **m);
//handling error
void nrerror(char error_text[]);
//complex operations
fcomplex cmul(fcomplex a, fcomplex b);
fcomplex cconj(fcomplex z);
fcomplex cadd(fcomplex a, fcomplex b);
float xcabs(fcomplex z);
float cphs(fcomplex z);
//functions
long next_pow2(long a);
void circ_shift(fcomplex *in, int na, int nc);
void left_shift(fcomplex *in, int na);
void right_shift(fcomplex *in, int na);
int roundfi(float a);
void sinc(int n, int m, float *coef);
void kaiser(int n, int m, float *coef, float beta);
void kaiser2(float beta, int n, float *coef);
void bandpass_filter(float bw, float bc, int n, int nfft, int ncshift, float beta, fcomplex *filter);
float bessi0(float x);
void four1(float data[], unsigned long nn, int isign);
//file operations
FILE *openfile(char *filename, char *pattern);
void readdata(void *data, size_t blocksize, FILE *fp);
void writedata(void *data, size_t blocksize, FILE *fp);
long file_length(FILE* fp, long width, long element_size);

View File

@ -0,0 +1,246 @@
//////////////////////////////////////
// Cunren Liang
// California Institute of Technology
// Copyright 2019
//////////////////////////////////////
//this program is tested against resamp.c, the outputs of the two are exactly the same.
#include "resamp.h"
//ALOS I or Q mean = 15.5, so get 15 or 16 randomly here
//#define ZERO_VALUE (char)(15 + rand() % 2)
//I changed the dynamic range when reading data
//ALOS I or Q mean = 63.5, so get 63 or 64 randomly here
#define ZERO_VALUE (char)(63 + rand() % 2)
typedef struct {
char re;
char im;
} char_complex;
char_complex *array1d_char_complex(long nc);
void free_array1d_char_complex(char_complex *fcv);
void normalize_kernel(float *kernel, long start_index, long end_index);
int resamp_azimuth(char *slc2, char *rslc2, int nrg, int naz1, int naz2, double prf, double *dopcoeff, double *azcoef, int n, double beta){
int i;
int verbose = 0;
if(verbose){
printf("\n\ninput parameters:\n");
printf("slc2: %s\n", slc2);
printf("rslc2: %s\n", rslc2);
printf("nrg: %d\n", nrg);
printf("naz1: %d\n", naz1);
printf("naz2: %d\n\n", naz2);
printf("prf: %f\n\n", prf);
for(i = 0; i < 4; i++){
printf("dopcoeff[%d]: %e\n", i, dopcoeff[i]);
}
printf("\n");
for(i = 0; i < 2; i++){
printf("azcoef[%d]: %e\n", i, azcoef[i]);
}
printf("\n");
}
FILE *slc2fp;
FILE *rslc2fp;
int m;
int interp_method;
int edge_method;
float azpos;
float azoff;
float az2;
int azi2;
float azf;
int azfn;
int hnm;
int hn;
float *sincc;
float *kaiserc;
float *kernel;
float *azkernel;
fcomplex *azkernel_fc;
fcomplex *rgrs;
fcomplex *azca;
fcomplex *rgrsb;
fcomplex *azrs;
char_complex *inl;
char_complex *outl;
float *dop;
float dopx;
fcomplex **inb;
int j, k, k1, k2;
int tmp1, tmp2;
int zero_flag;
float ftmp1, ftmp2;
fcomplex fctmp1, fctmp2;
m = 10000;
interp_method = 0;
edge_method = 2;
if((n % 2 == 0) || (n < 3)){
fprintf(stderr, "number of samples to be used in the resampling must be odd, and larger or equal to than 3\n");
exit(1);
}
slc2fp = openfile(slc2, "rb");
rslc2fp = openfile(rslc2, "wb");
hn = n / 2;
hnm = n * m / 2;
sincc = vector_float(-hnm, hnm);
kaiserc = vector_float(-hnm, hnm);
kernel = vector_float(-hnm, hnm);
azkernel = vector_float(-hn, hn);
azkernel_fc = vector_fcomplex(-hn, hn);
rgrs = vector_fcomplex(-hn, hn);
azca = vector_fcomplex(-hn, hn);
rgrsb = vector_fcomplex(-hn, hn);
azrs = array1d_fcomplex(nrg);
inl = array1d_char_complex(nrg);
outl = array1d_char_complex(nrg);
dop = array1d_float(nrg);
inb = array2d_fcomplex(naz2, nrg);
sinc(n, m, sincc);
kaiser(n, m, kaiserc, beta);
for(i = -hnm; i <= hnm; i++)
kernel[i] = kaiserc[i] * sincc[i];
for(i = 0; i < nrg; i++){
dop[i] = dopcoeff[0] + dopcoeff[1] * i + dopcoeff[2] * i * i + dopcoeff[3] * i * i * i;
if(verbose){
if(i % 500 == 0)
printf("range sample: %5d, doppler centroid frequency: %8.2f Hz\n", i, dop[i]);
}
}
for(i = 0; i < naz2; i++){
readdata((char_complex *)inl, (size_t)nrg * sizeof(char_complex), slc2fp);
for(j =0; j < nrg; j++){
inb[i][j].re = inl[j].re;
inb[i][j].im = inl[j].im;
}
}
for(i = 0; i < naz1; i++){
if((i + 1) % 100 == 0)
fprintf(stderr,"processing line: %6d of %6d\r", i+1, naz1);
for(j = 0; j < nrg; j++){
azrs[j].re = 0.0;
azrs[j].im = 0.0;
}
azpos = i;
azoff = azcoef[0] + azpos * azcoef[1];
az2 = i + azoff;
azi2 = roundfi(az2);
azf = az2 - azi2;
azfn = roundfi(azf * m);
if(edge_method == 0){
if(azi2 < hn || azi2 > naz2 - 1 - hn){
for(j = 0; j < nrg; j++){
outl[j].re = ZERO_VALUE;
outl[j].im = ZERO_VALUE;
}
writedata((char_complex *)outl, (size_t)nrg * sizeof(char_complex), rslc2fp);
continue;
}
}
else if(edge_method == 1){
if(azi2 < 0 || azi2 > naz2 - 1){
for(j = 0; j < nrg; j++){
outl[j].re = ZERO_VALUE;
outl[j].im = ZERO_VALUE;
}
writedata((char_complex *)outl, (size_t)nrg * sizeof(char_complex), rslc2fp);
continue;
}
}
else{
if(azi2 < -hn || azi2 > naz2 - 1 + hn){
for(j = 0; j < nrg; j++){
outl[j].re = ZERO_VALUE;
outl[j].im = ZERO_VALUE;
}
writedata((char_complex *)outl, (size_t)nrg * sizeof(char_complex), rslc2fp);
continue;
}
}
for(k = -hn; k <= hn; k++){
tmp2 = k * m - azfn;
if(tmp2 > hnm) tmp2 = hnm;
if(tmp2 < -hnm) tmp2 = -hnm;
azkernel[k] = kernel[tmp2];
}
normalize_kernel(azkernel, -hn, hn);
for(j = 0; j < nrg; j++){
for(k1 = -hn; k1 <= hn; k1++){
if((azi2 + k1 >= 0)&&(azi2 + k1 <= naz2-1)){
rgrs[k1].re = inb[azi2 + k1][j].re;
rgrs[k1].im = inb[azi2 + k1][j].im;
}
else{
rgrs[k1].re = ZERO_VALUE;
rgrs[k1].im = ZERO_VALUE;
}
}
dopx = dop[j];
for(k = -hn; k <= hn; k++){
ftmp1 = 2.0 * PI * dopx * k / prf;
azca[k].re = cos(ftmp1);
azca[k].im = sin(ftmp1);
if(interp_method == 0){
rgrsb[k] = cmul(rgrs[k], cconj(azca[k]));
azrs[j].re += rgrsb[k].re * azkernel[k];
azrs[j].im += rgrsb[k].im * azkernel[k];
}
else{
azkernel_fc[k].re = azca[k].re * azkernel[k];
azkernel_fc[k].im = azca[k].im * azkernel[k];
azrs[j] = cadd(azrs[j], cmul(rgrs[k], azkernel_fc[k]));
}
}
if(interp_method == 0){
ftmp1 = 2.0 * PI * dopx * azf / prf;
fctmp1.re = cos(ftmp1);
fctmp1.im = sin(ftmp1);
azrs[j] = cmul(azrs[j], fctmp1);
}
}
for(j = 0; j < nrg; j++){
outl[j].re = roundfi(azrs[j].re);
outl[j].im = roundfi(azrs[j].im);
}
writedata((char_complex *)outl, (size_t)nrg * sizeof(char_complex), rslc2fp);
}
fprintf(stderr,"processing line: %6d of %6d\n", naz1, naz1);
free_vector_float(sincc, -hnm, hnm);
free_vector_float(kaiserc, -hnm, hnm);
free_vector_float(kernel, -hnm, hnm);
free_vector_float(azkernel, -hn, hn);
free_vector_fcomplex(azkernel_fc, -hn, hn);
free_vector_fcomplex(rgrs, -hn, hn);
free_vector_fcomplex(azca, -hn, hn);
free_vector_fcomplex(rgrsb, -hn, hn);
free_array1d_fcomplex(azrs);
free_array1d_char_complex(inl);
free_array1d_char_complex(outl);
free_array1d_float(dop);
free_array2d_fcomplex(inb);
fclose(slc2fp);
fclose(rslc2fp);
return 0;
}
char_complex *array1d_char_complex(long nc){
char_complex *fcv;
fcv = (char_complex*) malloc(nc * sizeof(char_complex));
if(!fcv){
fprintf(stderr,"Error: cannot allocate 1-D char complex array\n");
exit(1);
}
return fcv;
}
void free_array1d_char_complex(char_complex *fcv){
free(fcv);
}
void normalize_kernel(float *kernel, long start_index, long end_index){
double sum;
long i;
sum = 0.0;
for(i = start_index; i <= end_index; i++)
sum += kernel[i];
if(sum!=0)
for(i = start_index; i <= end_index; i++)
kernel[i] /= sum;
}

View File

@ -1,48 +1,51 @@
//program updated to handle PRF change issue of ALOS-1
//Cunren Liang, December 2019
#include "image_sio.h"
#include "siocomplex.h"
#include <math.h>
fcomplex Cmul(fcomplex x, fcomplex y)
fcomplex_sio Cmul(fcomplex_sio x, fcomplex_sio y)
{
fcomplex z;
fcomplex_sio z;
z.r = x.r*y.r - x.i*y.i;
z.i = x.i*y.r + x.r*y.i;
return z;
}
fcomplex Cexp(float theta)
fcomplex_sio Cexp(float theta)
{
fcomplex z;
fcomplex_sio z;
z.r = cos(theta);
z.i = sin(theta);
return z;
}
fcomplex Conjg(fcomplex z)
fcomplex_sio Conjg(fcomplex_sio z)
{
fcomplex x;
fcomplex_sio x;
x.r = z.r;
x.i = -z.i;
return x;
}
fcomplex RCmul(float a, fcomplex z)
fcomplex_sio RCmul(float a, fcomplex_sio z)
{
fcomplex x;
fcomplex_sio x;
x.r = a*z.r;
x.i = a*z.i;
return x;
}
fcomplex Cadd(fcomplex x, fcomplex y)
fcomplex_sio Cadd(fcomplex_sio x, fcomplex_sio y)
{
fcomplex z;
fcomplex_sio z;
z.r = x.r + y.r;
z.i = x.i + y.i;
return z;
}
float Cabs(fcomplex z)
float Cabs(fcomplex_sio z)
{
return hypot(z.r, z.i);
}

View File

@ -1,11 +1,14 @@
//program updated to handle PRF change issue of ALOS-1
//Cunren Liang, December 2019
#ifndef _COMPLEX_H
#define _COMPLEX_H
fcomplex Cmul(fcomplex x, fcomplex y);
fcomplex Cexp(float theta);
fcomplex Conjg(fcomplex z);
fcomplex RCmul(float a, fcomplex z);
fcomplex Cadd(fcomplex x, fcomplex y);
float Cabs(fcomplex z);
fcomplex_sio Cmul(fcomplex_sio x, fcomplex_sio y);
fcomplex_sio Cexp(float theta);
fcomplex_sio Conjg(fcomplex_sio z);
fcomplex_sio RCmul(float a, fcomplex_sio z);
fcomplex_sio Cadd(fcomplex_sio x, fcomplex_sio y);
float Cabs(fcomplex_sio z);
#endif /* _COMPLEX_H */

View File

@ -25,7 +25,8 @@
// Author: Giangi Sacco
//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
//program updated to handle PRF change issue of ALOS-1
//Cunren Liang, December 2019
#include <Python.h>
@ -219,8 +220,8 @@ PyObject * ALOS_fbd2fbs_C(PyObject* self, PyObject* args)
i = j + r.first_sample;
/* increase dynamic range by 2 and set the mean value to 63.5 */
rtest = rintf(2.*cout[j].r+63.5);
itest = rintf(2.*cout[j].i+63.5);
rtest = rintf(cout[j].r+r.xmi);
itest = rintf(cout[j].i+r.xmq);
/* sometimes the range can exceed 0-127 so
clip the numbers to be in the correct range */

View File

@ -25,7 +25,8 @@
// Author: Giangi Sacco
//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
//program updated to handle PRF change issue of ALOS-1
//Cunren Liang, December 2019
#include <Python.h>
@ -197,10 +198,10 @@ PyObject * ALOS_fbs2fbd_C(PyObject* self, PyObject* args)
n4 = nffti/4;
for(i=0; i<n4;i++)
{
cout[i].r = cin[i].r;
cout[i].i = cin[i].i;
cout[i+n4].r = cin[i+3*n4].r;
cout[i+n4].i = cin[i+3*n4].i;
cout[i].r = 0.5*cin[i].r;
cout[i].i = 0.5*cin[i].i;
cout[i+n4].r = 0.5*cin[i+3*n4].r;
cout[i+n4].i = 0.5*cin[i+3*n4].i;
}
/*****Inverse FFT*****/
@ -219,8 +220,8 @@ PyObject * ALOS_fbs2fbd_C(PyObject* self, PyObject* args)
i = j + r.first_sample;
/* increase dynamic range by 2 and set the mean value to 63.5 */
rtest = rintf(2.*cout[j].r+63.5);
itest = rintf(2.*cout[j].i+63.5);
rtest = rintf(cout[j].r+r.xmi);
itest = rintf(cout[j].i+r.xmq);
/* sometimes the range can exceed 0-127 so
clip the numbers to be in the correct range */