[cig-commits] r18777 - mc/1D/hc/trunk
becker at geodynamics.org
becker at geodynamics.org
Sat Jul 16 12:19:08 PDT 2011
Author: becker
Date: 2011-07-16 12:19:07 -0700 (Sat, 16 Jul 2011)
New Revision: 18777
Added:
mc/1D/hc/trunk/hc_extract_spatial.c
Log:
*** empty log message ***
Added: mc/1D/hc/trunk/hc_extract_spatial.c
===================================================================
--- mc/1D/hc/trunk/hc_extract_spatial.c (rev 0)
+++ mc/1D/hc/trunk/hc_extract_spatial.c 2011-07-16 19:19:07 UTC (rev 18777)
@@ -0,0 +1,162 @@
+#include "hc.h"
+/*
+
+
+extract part of a spherical harmonics solution of a HC run
+
+
+$Id: hc_extract_sh_layer.c,v 1.9 2006/01/22 01:11:34 becker Exp becker $
+
+
+*/
+
+int main(int argc, char **argv)
+{
+ int ilayer,nsol,i,mode,shps,nset=1,loop,i1,i2;
+ FILE *in;
+ struct sh_lms *sol=NULL;
+ struct hcs *model;
+ HC_PREC fac[3] = {1.0,1.0,1.0};
+ hc_boolean binary = TRUE, verbose = FALSE, short_format = FALSE;
+ hc_struc_init(&model);
+ /*
+ deal with parameters
+ */
+ ilayer = 0;
+ mode = 1;
+ switch(argc){
+ case 3:
+ sscanf(argv[2],"%i",&ilayer);
+ break;
+ case 4:
+ sscanf(argv[2],"%i",&ilayer);
+ sscanf(argv[3],"%i",&mode);
+ break;
+ case 5:
+ sscanf(argv[2],"%i",&ilayer);
+ sscanf(argv[3],"%i",&mode);
+ sscanf(argv[4],"%i",&i);
+ if(i)
+ short_format = TRUE;
+ else
+ short_format = FALSE;
+ break;
+ default:
+ fprintf(stderr,"%s: usage\n%s sol.file layer [mode,%i] [short_format, %i]\n\n",
+ argv[0],argv[0],mode,short_format);
+ fprintf(stderr,"extracts spherical harmonic solution x (vel or str) from HC run\n");
+ fprintf(stderr,"layer: 1...nset\n");
+ fprintf(stderr,"\tif ilayer= 1..nset, will print one layer\n");
+ fprintf(stderr,"\t -1, will select nset\n");
+ fprintf(stderr,"\t -2, will print all layers\n");
+ fprintf(stderr,"mode: 1...6\n");
+ fprintf(stderr,"\tif mode = 1, will print x_r \n");
+ fprintf(stderr,"\t 2, will print x_pol x_tor \n");
+ fprintf(stderr,"\t 3, will print x_r x_pol x_tor\n");
+ fprintf(stderr,"\t 4, will print the depth levels of all layers\n");
+ fprintf(stderr,"\t 5, will print x_pol\n");
+ fprintf(stderr,"\t 6, will print x_tor\n");
+ exit(-1);
+ break;
+ }
+ if(mode == 4)
+ ilayer = -2;
+ /*
+ read in solution
+ */
+ in = ggrd_open(argv[1],"r","hc_extract_sh_layer");
+ hc_read_sh_solution(model,&sol,in,binary,verbose);
+ fclose(in);
+ nsol = model->nradp2 * 3;
+ /*
+ deal with selection
+ */
+ loop = 0;
+ if(ilayer == -1)
+ ilayer = model->nradp2;
+ if(ilayer == -2){
+ ilayer = model->nradp2;
+ loop =1;
+ }
+ if((ilayer<1)||(ilayer > model->nradp2)){
+ fprintf(stderr,"%s: ilayer (%i) out of range, use 1 ... %i\n",
+ argv[0],ilayer,model->nradp2);
+ exit(-1);
+ }
+ if(loop){
+ i1=0;i2=model->nradp2-1;
+ if(short_format)
+ printf("%i\n",model->nradp2);
+ }else{
+ i1=ilayer-1;i2 = i1;
+ }
+ /* detect number of expansions */
+ if((mode == 1)||(mode == 5)||(mode == 6))
+ shps = 1;
+ else if(mode == 2)
+ shps = 2;
+ else if(mode == 3)
+ shps = 3;
+ for(ilayer=i1;ilayer <= i2;ilayer++){
+ /*
+ output
+ */
+ if(mode != 4){
+ /* SH header */
+ if(short_format && loop)
+ fprintf(stdout,"%g\n",HC_Z_DEPTH(model->r[ilayer]));
+ sh_print_parameters_to_file((sol+ilayer*3),shps,
+ ilayer,nset,(HC_PREC)(HC_Z_DEPTH(model->r[ilayer])),
+ stdout,short_format,FALSE,verbose);
+ }
+ switch(mode){
+ case 1:
+ /* */
+ if(verbose)
+ fprintf(stderr,"%s: printing x_r SHE at layer %i (depth: %g)\n",
+ argv[0],ilayer,HC_Z_DEPTH(model->r[ilayer]));
+ sh_print_coefficients_to_file((sol+ilayer*3),shps,stdout,fac,FALSE,verbose);
+ break;
+ case 2:
+ /* */
+ if(verbose)
+ fprintf(stderr,"%s: printing x_pol x_tor SHE at layer %i (depth: %g)\n",
+ argv[0],ilayer,HC_Z_DEPTH(model->r[ilayer]));
+ sh_print_coefficients_to_file((sol+ilayer*3+1),shps,stdout,fac,FALSE,verbose);
+ break;
+ case 3:
+ /* mode == 3 */
+ if(verbose)
+ fprintf(stderr,"%s: printing x_r x_pol x_tor SHE at layer %i (depth: %g)\n",
+ argv[0],ilayer,HC_Z_DEPTH(model->r[ilayer]));
+ sh_print_coefficients_to_file((sol+ilayer*3),shps,stdout,fac,FALSE,verbose);
+ break;
+ case 4:
+ fprintf(stdout,"%5i %11g\n",ilayer,HC_Z_DEPTH(model->r[ilayer]));
+ break;
+ case 5:
+ /* */
+ if(verbose)
+ fprintf(stderr,"%s: printing x_pol SHE at layer %i (depth: %g)\n",
+ argv[0],ilayer,HC_Z_DEPTH(model->r[ilayer]));
+ sh_print_coefficients_to_file((sol+ilayer*3+1),shps,stdout,fac,FALSE,verbose);
+ break;
+ case 6:
+ /* */
+ if(verbose)
+ fprintf(stderr,"%s: printing x_tor SHE at layer %i (depth: %g)\n",
+ argv[0],ilayer,HC_Z_DEPTH(model->r[ilayer]));
+ sh_print_coefficients_to_file((sol+ilayer*3+2),shps,stdout,fac,FALSE,verbose);
+ break;
+
+ default:
+ fprintf(stderr,"%s: error, mode %i undefined\n",argv[0],mode);
+ exit(-1);
+ break;
+ }
+ }
+ /* clear and exit */
+ sh_free_expansion(sol,nsol);
+
+ return 0;
+}
More information about the CIG-COMMITS
mailing list