Commit b0f98aaf authored by Florian Wittkamp's avatar Florian Wittkamp

Removed zero-phase option from the code

I removed the zero phase option from the code, due to complications depending on the source signal and the time window. However, from a technically point of view it was implemented right.
parent 42299882
...@@ -169,7 +169,6 @@ ...@@ -169,7 +169,6 @@
"F_LOW_PASS_END" : "75.0", "F_LOW_PASS_END" : "75.0",
"F_LOW_PASS_INCR" : "10.0", "F_LOW_PASS_INCR" : "10.0",
"ORDER" : "2", "ORDER" : "2",
"ZERO_PHASE" : "0",
"FREQ_FILE" : "frequencies.dat", "FREQ_FILE" : "frequencies.dat",
"WRITE_FILTERED_DATA" : "0", "WRITE_FILTERED_DATA" : "0",
......
...@@ -59,7 +59,7 @@ void exchange_par(void){ ...@@ -59,7 +59,7 @@ void exchange_par(void){
extern float npower, k_max_PML; extern float npower, k_max_PML;
extern int INV_STF, N_STF, N_STF_START; extern int INV_STF, N_STF, N_STF_START;
extern char PARA[STRING_SIZE]; extern char PARA[STRING_SIZE];
extern int TIME_FILT, ORDER, ZERO_PHASE,WRITE_FILTERED_DATA; extern int TIME_FILT, ORDER,WRITE_FILTERED_DATA;
extern float F_LOW_PASS_START, F_LOW_PASS_END, F_LOW_PASS_INCR, F_HIGH_PASS; extern float F_LOW_PASS_START, F_LOW_PASS_END, F_LOW_PASS_INCR, F_HIGH_PASS;
extern int LNORM, DTINV; extern int LNORM, DTINV;
extern int STEPMAX; extern int STEPMAX;
...@@ -325,7 +325,7 @@ void exchange_par(void){ ...@@ -325,7 +325,7 @@ void exchange_par(void){
idum[85] = NO_OF_TESTSHOTS; idum[85] = NO_OF_TESTSHOTS;
idum[86] = ZERO_PHASE; // idum[86] = EMPTY;
idum[87] = VELOCITY; idum[87] = VELOCITY;
...@@ -609,7 +609,7 @@ void exchange_par(void){ ...@@ -609,7 +609,7 @@ void exchange_par(void){
NO_OF_TESTSHOTS = idum[85]; NO_OF_TESTSHOTS = idum[85];
ZERO_PHASE = idum[86]; // EMPTY = idum[86];
VELOCITY = idum[87]; VELOCITY = idum[87];
......
...@@ -72,7 +72,7 @@ int TAPER_STF; ...@@ -72,7 +72,7 @@ int TAPER_STF;
int INV_STF, N_STF, N_STF_START,TW_IND; int INV_STF, N_STF, N_STF_START,TW_IND;
char PARA[STRING_SIZE]; char PARA[STRING_SIZE];
int TIME_FILT, ORDER, ZERO_PHASE; int TIME_FILT, ORDER;
int WRITE_FILTERED_DATA; int WRITE_FILTERED_DATA;
float F_LOW_PASS_START, F_LOW_PASS_END, F_LOW_PASS_INCR, F_HIGH_PASS; float F_LOW_PASS_START, F_LOW_PASS_END, F_LOW_PASS_INCR, F_HIGH_PASS;
......
...@@ -62,7 +62,7 @@ void read_par_json(FILE *fp, char *fileinp){ ...@@ -62,7 +62,7 @@ void read_par_json(FILE *fp, char *fileinp){
extern int INV_STF, N_STF, N_STF_START; extern int INV_STF, N_STF, N_STF_START;
extern char PARA[STRING_SIZE]; extern char PARA[STRING_SIZE];
extern int TIME_FILT, ORDER, ZERO_PHASE,WRITE_FILTERED_DATA; extern int TIME_FILT, ORDER,WRITE_FILTERED_DATA;
extern float F_LOW_PASS_START, F_LOW_PASS_END, F_LOW_PASS_INCR, F_HIGH_PASS; extern float F_LOW_PASS_START, F_LOW_PASS_END, F_LOW_PASS_INCR, F_HIGH_PASS;
extern int LNORM, DTINV; extern int LNORM, DTINV;
...@@ -847,9 +847,6 @@ void read_par_json(FILE *fp, char *fileinp){ ...@@ -847,9 +847,6 @@ void read_par_json(FILE *fp, char *fileinp){
declare_error("Variable ORDER could not be retrieved from the json input file!"); declare_error("Variable ORDER could not be retrieved from the json input file!");
} }
} }
if (get_int_from_objectlist("ZERO_PHASE",number_readobjects,&ZERO_PHASE,varname_list, value_list)){
ZERO_PHASE=0;
fprintf(fp,"Variable ZERO_PHASE is set to default value %i.\n",ZERO_PHASE);}
if (TIME_FILT==2) { if (TIME_FILT==2) {
if (get_float_from_objectlist("F_HIGH_PASS",number_readobjects,&F_HIGH_PASS,varname_list, value_list)){ if (get_float_from_objectlist("F_HIGH_PASS",number_readobjects,&F_HIGH_PASS,varname_list, value_list)){
F_HIGH_PASS=0.0; F_HIGH_PASS=0.0;
...@@ -858,9 +855,6 @@ void read_par_json(FILE *fp, char *fileinp){ ...@@ -858,9 +855,6 @@ void read_par_json(FILE *fp, char *fileinp){
declare_error("Variable FREQ_FILE could not be retrieved from the json input file!"); declare_error("Variable FREQ_FILE could not be retrieved from the json input file!");
if (get_int_from_objectlist("ORDER",number_readobjects,&ORDER,varname_list, value_list)) if (get_int_from_objectlist("ORDER",number_readobjects,&ORDER,varname_list, value_list))
declare_error("Variable ORDER could not be retrieved from the json input file!"); declare_error("Variable ORDER could not be retrieved from the json input file!");
if (get_int_from_objectlist("ZERO_PHASE",number_readobjects,&ZERO_PHASE,varname_list, value_list)){
ZERO_PHASE=0;
fprintf(fp,"Variable ZERO_PHASE is set to default value %i.\n",ZERO_PHASE);}
} }
} }
......
...@@ -39,7 +39,7 @@ void timedomain_filt(float ** data, float fc, int order, int ntr, int ns, int m ...@@ -39,7 +39,7 @@ void timedomain_filt(float ** data, float fc, int order, int ntr, int ns, int m
/* declaration of extern variables */ /* declaration of extern variables */
extern float DT, F_HIGH_PASS; extern float DT, F_HIGH_PASS;
extern int ZERO_PHASE, NT,MYID; extern int NT,MYID;
/* declaration of local variables */ /* declaration of local variables */
int itr, j, ns_reverse; int itr, j, ns_reverse;
...@@ -47,10 +47,8 @@ void timedomain_filt(float ** data, float fc, int order, int ntr, int ns, int m ...@@ -47,10 +47,8 @@ void timedomain_filt(float ** data, float fc, int order, int ntr, int ns, int m
double *seismogram_hp, *seismogram_reverse_hp, T0_hp; double *seismogram_hp, *seismogram_reverse_hp, T0_hp;
seismogram = dvector(1,ns); seismogram = dvector(1,ns);
if (ZERO_PHASE==1) seismogram_reverse = dvector(1,ns);
seismogram_hp = dvector(1,ns); seismogram_hp = dvector(1,ns);
if (ZERO_PHASE==1) seismogram_reverse_hp = dvector(1,ns);
T0=1.0/(double)fc; T0=1.0/(double)fc;
if(F_HIGH_PASS) if(F_HIGH_PASS)
...@@ -65,18 +63,6 @@ void timedomain_filt(float ** data, float fc, int order, int ntr, int ns, int m ...@@ -65,18 +63,6 @@ void timedomain_filt(float ** data, float fc, int order, int ntr, int ns, int m
} }
seife_lpb(seismogram,ns+1,DT,T0,order); /* ns+1 because vector[0] is also allocated and otherwise seife_lpb do not filter the last sample */ seife_lpb(seismogram,ns+1,DT,T0,order); /* ns+1 because vector[0] is also allocated and otherwise seife_lpb do not filter the last sample */
if (ZERO_PHASE==1){
ns_reverse=ns;
for (j=1;j<=ns;j++) {
seismogram_reverse[ns_reverse]=seismogram[j];
ns_reverse--;}
seife_lpb(seismogram_reverse,ns+1,DT,T0,order);
ns_reverse=ns;
for (j=1;j<=ns;j++) {
seismogram[ns_reverse]=seismogram_reverse[j];
ns_reverse--;}
}
for (j=1;j<=ns;j++){ for (j=1;j<=ns;j++){
data[itr][j]=(float)seismogram[j]; data[itr][j]=(float)seismogram[j];
...@@ -92,17 +78,6 @@ void timedomain_filt(float ** data, float fc, int order, int ntr, int ns, int m ...@@ -92,17 +78,6 @@ void timedomain_filt(float ** data, float fc, int order, int ntr, int ns, int m
seife_hpb(seismogram_hp,ns+1,DT,T0_hp,order); seife_hpb(seismogram_hp,ns+1,DT,T0_hp,order);
if (ZERO_PHASE==1){
ns_reverse=ns;
for (j=1;j<=ns;j++) {
seismogram_reverse_hp[ns_reverse]=seismogram_hp[j];
ns_reverse--;}
seife_hpb(seismogram_reverse_hp,ns+1,DT,T0_hp,order);
ns_reverse=ns;
for (j=1;j<=ns;j++) {
seismogram_hp[ns_reverse]=seismogram_reverse_hp[j];
ns_reverse--;}
}
for (j=1;j<=ns;j++){ for (j=1;j<=ns;j++){
data[itr][j]=(float)seismogram_hp[j]; data[itr][j]=(float)seismogram_hp[j];
} }
...@@ -110,9 +85,7 @@ void timedomain_filt(float ** data, float fc, int order, int ntr, int ns, int m ...@@ -110,9 +85,7 @@ void timedomain_filt(float ** data, float fc, int order, int ntr, int ns, int m
} /* end of itr<=ntr loop */ } /* end of itr<=ntr loop */
free_dvector(seismogram,1,ns); free_dvector(seismogram,1,ns);
if (ZERO_PHASE==1) free_dvector(seismogram_reverse,1,ns);
free_dvector(seismogram_hp,1,ns); free_dvector(seismogram_hp,1,ns);
if (ZERO_PHASE==1) free_dvector(seismogram_reverse_hp,1,ns);
} /* end of function */ } /* end of function */
...@@ -39,7 +39,7 @@ void timedomain_filt_vector(float * data, float fc, int order, int ns, int meth ...@@ -39,7 +39,7 @@ void timedomain_filt_vector(float * data, float fc, int order, int ns, int meth
/* declaration of external variables */ /* declaration of external variables */
extern float DT, F_HIGH_PASS; extern float DT, F_HIGH_PASS;
extern int ZERO_PHASE, NT,MYID; extern int NT,MYID;
/* declaration of local variables */ /* declaration of local variables */
int itr, j, ns_reverse; int itr, j, ns_reverse;
...@@ -50,10 +50,8 @@ void timedomain_filt_vector(float * data, float fc, int order, int ns, int meth ...@@ -50,10 +50,8 @@ void timedomain_filt_vector(float * data, float fc, int order, int ns, int meth
declare_error("Order of timedomain filter must be an even number!");} declare_error("Order of timedomain filter must be an even number!");}
seismogram = dvector(1,ns); seismogram = dvector(1,ns);
if (ZERO_PHASE==1) seismogram_reverse = dvector(1,ns);
seismogram_hp = dvector(1,ns); seismogram_hp = dvector(1,ns);
if (ZERO_PHASE==1) seismogram_reverse_hp = dvector(1,ns);
T0=1.0/(double)fc; T0=1.0/(double)fc;
if(F_HIGH_PASS) if(F_HIGH_PASS)
...@@ -67,18 +65,6 @@ void timedomain_filt_vector(float * data, float fc, int order, int ns, int meth ...@@ -67,18 +65,6 @@ void timedomain_filt_vector(float * data, float fc, int order, int ns, int meth
} }
seife_lpb(seismogram,ns+1,DT,T0,order); /* ns+1 because vector[0] is also allocated and otherwise seife_lpb do not filter the last sample */ seife_lpb(seismogram,ns+1,DT,T0,order); /* ns+1 because vector[0] is also allocated and otherwise seife_lpb do not filter the last sample */
if (ZERO_PHASE==1){
ns_reverse=ns;
for (j=1;j<=ns;j++) {
seismogram_reverse[ns_reverse]=seismogram[j];
ns_reverse--;}
seife_lpb(seismogram_reverse,ns+1,DT,T0,order);
ns_reverse=ns;
for (j=1;j<=ns;j++) {
seismogram[ns_reverse]=seismogram_reverse[j];
ns_reverse--;}
}
for (j=1;j<=ns;j++){ for (j=1;j<=ns;j++){
data[j]=(float)seismogram[j]; data[j]=(float)seismogram[j];
...@@ -92,26 +78,13 @@ void timedomain_filt_vector(float * data, float fc, int order, int ns, int meth ...@@ -92,26 +78,13 @@ void timedomain_filt_vector(float * data, float fc, int order, int ns, int meth
seife_hpb(seismogram_hp,ns+1,DT,T0_hp,order); seife_hpb(seismogram_hp,ns+1,DT,T0_hp,order);
if (ZERO_PHASE==1){
ns_reverse=ns;
for (j=1;j<=ns;j++) {
seismogram_reverse_hp[ns_reverse]=seismogram_hp[j];
ns_reverse--;}
seife_hpb(seismogram_reverse_hp,ns+1,DT,T0_hp,order);
ns_reverse=ns;
for (j=1;j<=ns;j++) {
seismogram_hp[ns_reverse]=seismogram_reverse_hp[j];
ns_reverse--;}
}
for (j=1;j<=ns;j++){ for (j=1;j<=ns;j++){
data[j]=(float)seismogram_hp[j]; data[j]=(float)seismogram_hp[j];
} }
} }
free_dvector(seismogram,1,ns); free_dvector(seismogram,1,ns);
if (ZERO_PHASE==1) free_dvector(seismogram_reverse,1,ns);
free_dvector(seismogram_hp,1,ns); free_dvector(seismogram_hp,1,ns);
if (ZERO_PHASE==1) free_dvector(seismogram_reverse_hp,1,ns);
} }
\ No newline at end of file
...@@ -62,7 +62,7 @@ void write_par(FILE *fp){ ...@@ -62,7 +62,7 @@ void write_par(FILE *fp){
extern int INV_STF, N_STF, N_STF_START; extern int INV_STF, N_STF, N_STF_START;
extern char PARA[STRING_SIZE]; extern char PARA[STRING_SIZE];
extern int TIME_FILT, ORDER, ZERO_PHASE; extern int TIME_FILT, ORDER;
extern float F_LOW_PASS_START, F_LOW_PASS_END, F_LOW_PASS_INCR; extern float F_LOW_PASS_START, F_LOW_PASS_END, F_LOW_PASS_INCR;
extern int LNORM, DTINV; extern int LNORM, DTINV;
extern int STEPMAX, TRKILL, TRKILL_STF; extern int STEPMAX, TRKILL, TRKILL_STF;
...@@ -525,11 +525,9 @@ void write_par(FILE *fp){ ...@@ -525,11 +525,9 @@ void write_par(FILE *fp){
fprintf(fp," Order of lowpass filter is:\t%d\n",ORDER); fprintf(fp," Order of lowpass filter is:\t%d\n",ORDER);
if ((ORDER%2)!=0){ if ((ORDER%2)!=0){
declare_error(" Order of time domain filter must be an even number! \n");} declare_error(" Order of time domain filter must be an even number! \n");}
if (ZERO_PHASE){ } else {
fprintf(fp, " ZERO_PHASE=%i: Zero phase filtering is applied! \n",ZERO_PHASE);} fprintf(fp," TIME_FILT=%d: No time domain filtering is applied.\n",TIME_FILT);
else fprintf(fp, " No zero phase filtering is applied! \n");} }
else fprintf(fp," TIME_FILT=%d: No time domain filtering is applied.\n",TIME_FILT);
fprintf(fp,"\n\n"); fprintf(fp,"\n\n");
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment