Changeset 5477


Ignore:
Timestamp:
Jul 7, 2008, 6:02:36 PM (11 years ago)
Author:
jakeman
Message:

cleaned up urs_ext.c. Currently failing to read fros and lros correctly in in boxing_day mux2 files

File:
1 edited

Legend:

Unmodified
Added
Removed
  • anuga_core/source/anuga/shallow_water/urs_ext.c

    r5476 r5477  
    144144   int len_sts_data;
    145145   float *temp_sts_data;
     146
    146147   PyArrayObject *sts_data;
    147148   
     
    149150   PyArrayObject *permutation_array;
    150151   
    151    long int offset;
    152 
    153    /* Allocate space for the names and the weights and pointers to the data*/   
     152   long int offset;   
    154153   
    155154   /* Check that the input files have mux2 extension*/
     
    187186   fread(&nsta0,sizeof(int),1,fp);
    188187   
    189    if (permutation == Py_None){
    190      // if no permutation is specified return the times series of all the gauges in the mux file
    191      permuation_dimensions[0]=nsta0;
    192      permutation_array=(PyArrayObject *) PyArray_FromDims(1, permuation_dimensions, PyArray_INT);
    193      for (ista=0; ista<nsta0; ista++){
    194        *(long *) (permutation_array -> data+ista*permutation_array->strides[0])=ista;
    195      }
    196    }else{
    197      // Specifies the gauge numbers that for which data is to be extracted
    198      permutation_array=(PyArrayObject *) PyArray_ContiguousFromObject(permutation,PyArray_INT,1,1);
    199      N = permutation_array->dimensions[0]; //FIXME: this is overwritten below
    200      for (ista=0; ista<N; ista++){
    201        if ((int)*(long *) (permutation_array -> data+ista*permutation_array->strides[0])>=nsta0){
    202          printf("Maximum index = %d, you had %d\n",nsta0-1,(int)*(long *) (permutation_array -> data+ista*permutation_array->strides[0]));
    203          PyErr_SetString(PyExc_RuntimeError,"The permutation specified is out of bounds");
    204          return NULL;
    205        }
    206      }
    207    }
    208    if(permutation_array == NULL){
    209      PyErr_SetString(PyExc_RuntimeError,"Memory for permutation_array array could not be allocated");
    210      return NULL;
    211    }
    212188   /*now allocate space for, and read in, the structures for each station*/
    213189   mytgs0 = (struct tgsrwg *) malloc(nsta0*sizeof(struct tgsrwg));
     
    216192   /*make an array to hold the start and stop steps for each station for each
    217193     source*/   
    218    //FIXME: Should only store start and stop times for one source
    219194   fros = (int *) malloc(nsta0*numSrc*sizeof(int));
    220195   lros = (int *) malloc(nsta0*numSrc*sizeof(int));
     
    235210   }
    236211   fclose(fp);
     212
    237213   // Allocate header array for each remaining source file.
    238    // FIXME: only need to store for one source and which is compared to all subsequent
    239    // source headers
    240214   if(numSrc > 1){
    241215      /* allocate space for tgsrwg for the other sources */
     
    243217   } else {
    244218     /* FIXME (JJ): What should happen in case the are no source files?*/
    245      /* If we exit here, tests will break */
    246      // fprintf(stderr, "No source file specified\n");
    247      // exit(-1);       
    248    }
    249    
    250    /* loop over remaining sources, check compatibility, and read them into
    251     *muxData */
     219     /* If we exit here, tests will break */       
     220   }
     221   
     222   /* loop over remaining sources and read headers */
    252223   for(isrc=1; isrc<numSrc; isrc++){
    253224     muxFileName = muxFileNameArray[isrc];
     
    291262      numData = getNumData(fros+isrc*nsta0, lros+isrc*nsta0, nsta0);
    292263
    293       /* Burbidge: Added a sanity check here */
    294264      if (numData < 0){
    295265          fprintf(stderr,"Size of data block appears to be negative!\n");
     
    299269      fclose(fp);             
    300270   }
     271
     272   for (i=0;i<nsta0*numSrc;i++){
     273       printf("%d, fros %d\n",i,fros[i]);
     274       printf("%d, lros %d\n",i,lros[i]);
     275   }
     276
     277   if (permutation == Py_None){
     278       // if no permutation is specified return the times series of all the gauges in the mux file
     279       permuation_dimensions[0]=nsta0;
     280       permutation_array=(PyArrayObject *) PyArray_FromDims(1, permuation_dimensions, PyArray_INT);
     281       for (ista=0; ista<nsta0; ista++){
     282           *(long *) (permutation_array -> data+ista*permutation_array->strides[0])=ista;
     283       }
     284   }else{
     285       // Specifies the gauge numbers that for which data is to be extracted
     286       permutation_array=(PyArrayObject *) PyArray_ContiguousFromObject(permutation,PyArray_INT,1,1);
     287       N = permutation_array->dimensions[0]; //FIXME: this is overwritten below
     288       for (ista=0; ista<N; ista++){
     289           if ((int)*(long *) (permutation_array -> data+ista*permutation_array->strides[0])>=nsta0){
     290               printf("Maximum index = %d, you had %d\n",nsta0-1,(int)*(long *) (permutation_array -> data+ista*permutation_array->strides[0]));
     291               PyErr_SetString(PyExc_RuntimeError,"The permutation specified is out of bounds");
     292               return NULL;
     293           }
     294       }
     295   }
     296   if(permutation_array == NULL){
     297     PyErr_SetString(PyExc_RuntimeError,"Memory for permutation_array array could not be allocated");
     298     return NULL;
     299   }
     300
    301301   N = permutation_array->dimensions[0];
    302302   
     
    307307   
    308308   // Find min and max start times of all gauges
    309    // FIXME: ONLY CHECK gauges in permutation for start and finsh times
    310309   start_tstep=num_ts+1;
    311310   finish_tstep=-1;
    312    for (ista=0;ista<nsta0;ista++){
    313      if (fros[ista]< start_tstep){
    314        start_tstep=fros[ista];
    315      }
    316      if (lros[0+nsta0*ista] > finish_tstep){
    317        finish_tstep=lros[ista];
    318      }
     311   for (i=0;i<N;i++){
     312       ista = *(long *) (permutation_array -> data+i*permutation_array->strides[0]); 
     313       if (fros[ista]< start_tstep){
     314           start_tstep=fros[ista];
     315       }
     316       if (lros[0+nsta0*ista] > finish_tstep){
     317           finish_tstep=lros[ista];
     318       }
    319319   }
    320320   
    321321   if ((start_tstep>num_ts) | (finish_tstep < 0)){
     322       printf("s=%d,f=%d,num_ts=%d\n",start_tstep,finish_tstep,num_ts);
    322323     PyErr_SetString(PyExc_RuntimeError,"Gauge data has incorrect start and finsh times");
    323324     return NULL;
     
    482483}
    483484
    484 /* Burbidge: No "offset" is sent. Replace with max. Added grid_id */
    485 void wrttxt( fname, dt, nt, x, beg, lat, lon, max, depth, grid_id, istart, istop )
    486 char *fname;
    487 float dt, *x, beg, max, lat, lon, depth;
    488 int grid_id;
    489 int nt;
    490 int istart, istop;
    491 {
    492    int it;
    493    float t;
    494    FILE *fp;
    495 
    496    fp = fopen(fname,"w");
    497    fprintf(fp,"> lat %g lon %g depth %g max %g start_time %g stop_time %g grid_id %d\n", lat, lon, depth, max,
    498             istart*dt, istop*dt, grid_id );
    499    for(it=0;it<nt;it++)
    500    {
    501       t=beg+it*dt;
    502       fprintf(fp,"%9.3f %g\n", t, x[it]);
    503    }
    504    fclose(fp);
    505 }
    506  
     485 
    507486char isdata(float x)
    508487{
Note: See TracChangeset for help on using the changeset viewer.