Changeset 5477
 Timestamp:
 Jul 7, 2008, 6:02:36 PM (11 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

anuga_core/source/anuga/shallow_water/urs_ext.c
r5476 r5477 144 144 int len_sts_data; 145 145 float *temp_sts_data; 146 146 147 PyArrayObject *sts_data; 147 148 … … 149 150 PyArrayObject *permutation_array; 150 151 151 long int offset; 152 153 /* Allocate space for the names and the weights and pointers to the data*/ 152 long int offset; 154 153 155 154 /* Check that the input files have mux2 extension*/ … … 187 186 fread(&nsta0,sizeof(int),1,fp); 188 187 189 if (permutation == Py_None){190 // if no permutation is specified return the times series of all the gauges in the mux file191 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 extracted198 permutation_array=(PyArrayObject *) PyArray_ContiguousFromObject(permutation,PyArray_INT,1,1);199 N = permutation_array>dimensions[0]; //FIXME: this is overwritten below200 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",nsta01,(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 }212 188 /*now allocate space for, and read in, the structures for each station*/ 213 189 mytgs0 = (struct tgsrwg *) malloc(nsta0*sizeof(struct tgsrwg)); … … 216 192 /*make an array to hold the start and stop steps for each station for each 217 193 source*/ 218 //FIXME: Should only store start and stop times for one source219 194 fros = (int *) malloc(nsta0*numSrc*sizeof(int)); 220 195 lros = (int *) malloc(nsta0*numSrc*sizeof(int)); … … 235 210 } 236 211 fclose(fp); 212 237 213 // Allocate header array for each remaining source file. 238 // FIXME: only need to store for one source and which is compared to all subsequent239 // source headers240 214 if(numSrc > 1){ 241 215 /* allocate space for tgsrwg for the other sources */ … … 243 217 } else { 244 218 /* 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 */ 252 223 for(isrc=1; isrc<numSrc; isrc++){ 253 224 muxFileName = muxFileNameArray[isrc]; … … 291 262 numData = getNumData(fros+isrc*nsta0, lros+isrc*nsta0, nsta0); 292 263 293 /* Burbidge: Added a sanity check here */294 264 if (numData < 0){ 295 265 fprintf(stderr,"Size of data block appears to be negative!\n"); … … 299 269 fclose(fp); 300 270 } 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",nsta01,(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 301 301 N = permutation_array>dimensions[0]; 302 302 … … 307 307 308 308 // Find min and max start times of all gauges 309 // FIXME: ONLY CHECK gauges in permutation for start and finsh times310 309 start_tstep=num_ts+1; 311 310 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 } 319 319 } 320 320 321 321 if ((start_tstep>num_ts)  (finish_tstep < 0)){ 322 printf("s=%d,f=%d,num_ts=%d\n",start_tstep,finish_tstep,num_ts); 322 323 PyErr_SetString(PyExc_RuntimeError,"Gauge data has incorrect start and finsh times"); 323 324 return NULL; … … 482 483 } 483 484 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 507 486 char isdata(float x) 508 487 {
Note: See TracChangeset
for help on using the changeset viewer.