Changeset 5612


Ignore:
Timestamp:
Aug 5, 2008, 3:46:33 PM (15 years ago)
Author:
ole
Message:

A bit more cleanup in urs2sts - all tests pass, but I think we need a few more.

Location:
anuga_core/source/anuga/shallow_water
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • anuga_core/source/anuga/shallow_water/data_manager.py

    r5611 r5612  
    48884888    numSrc=len(filenames)
    48894889
    4890     file_params=-1*ones(3,Float)#[nsta,dt,nt]
     4890    file_params=-1*ones(3,Float) #[nsta,dt,nt]
    48914891   
    48924892    # Convert verbose to int C flag
    … …  
    48964896        verbose=0
    48974897       
    4898        
    4899     #msg = 'I got the permutation vector:' + str(permutation)
    4900     #assert permutation is not None, msg
    49014898    if permutation is None:
    49024899        permutation = ensure_numeric([], Float)   
    49034900
    49044901    # Call underlying C implementation urs2sts_ext.c   
    4905     data = read_mux2(numSrc,filenames,weights,file_params,permutation,verbose)
     4902    data = read_mux2(numSrc, filenames, weights, file_params, permutation, verbose)
    49064903
    49074904    msg='File parameter values were not read in correctly from c file'
    … …  
    49254922    msg='Must have at least one gauge value'
    49264923    assert nt>0,msg
    4927     
     4924   
    49284925    OFFSET=5 # Number of site parameters p passed back with data
    49294926             # p=[geolat,geolon,depth,start_tstep,finish_tstep]
    4930 
     4927     
     4928    # FIXME (Ole): What is the relationship with params and data.shape ?
     4929    # It looks as if the following asserts should pass but they don't always
     4930    #
     4931             
     4932    #print
     4933    #print nsta, nt, data.shape
     4934       
     4935    #msg = 'nt = %d, data.shape[1] == %d' %(nt, data.shape[1])
     4936    #assert nt == data.shape[1] - OFFSET, msg
     4937   
     4938    #msg = 'nsta = %d, data.shape[0] == %d' %(nsta, data.shape[0])   
     4939    #assert nsta == data.shape[0], msg
     4940
     4941   
     4942    # Number of stations in ordering file
     4943    number_of_selected_stations = data.shape[0]
     4944
     4945    # Index where data ends and parameters begin
    49314946    parameters_index = data.shape[1]-OFFSET         
    49324947             
    49334948    times=dt*arange(parameters_index)   
    4934     latitudes=zeros(data.shape[0], Float)
    4935     longitudes=zeros(data.shape[0], Float)
    4936     elevation=zeros(data.shape[0], Float)
    4937     quantity=zeros((data.shape[0], parameters_index), Float)
     4949    latitudes=zeros(number_of_selected_stations, Float)
     4950    longitudes=zeros(number_of_selected_stations, Float)
     4951    elevation=zeros(number_of_selected_stations, Float)
     4952    quantity=zeros((number_of_selected_stations, parameters_index), Float)
    49384953   
    49394954   
    49404955    starttime=1e16
    4941     for i in range(0, data.shape[0]):
     4956    for i in range(number_of_selected_stations):
     4957        quantity[i][:]=data[i][:parameters_index]
     4958   
    49424959        latitudes[i]=data[i][parameters_index]
    49434960        longitudes[i]=data[i][parameters_index+1]
    49444961        elevation[i]=-data[i][parameters_index+2]
    4945         quantity[i][:]=data[i][:parameters_index] # Was data[i][:-OFFSET]
    4946        
    49474962        first_time_step = data[i][parameters_index+3]
    4948         #print 'datamanager:', i, first_time_step, dt*first_time_step
    4949         starttime=min(dt*data[i][parameters_index+3],starttime)
     4963       
     4964        starttime=min(dt*first_time_step, starttime)
    49504965       
    49514966    return times, latitudes, longitudes, elevation, quantity, starttime

  • anuga_core/source/anuga/shallow_water/urs_ext.c

    r5611 r5612  
    203203    }
    204204
    205     /* Loop over all sources, read headers and check compatibility */
     205    // Loop over all sources, read headers and check compatibility
    206206    for (i = 0; i < numSrc; i++)
    207207    {
    208208        muxFileName = muxFileNameArray[i];
    209209
    210         /* open the mux file */
     210        // Open the mux file
    211211        if((fp = fopen(muxFileName, "r")) == NULL)
    212212        {
    … …  
    219219            fread(total_number_of_stations, sizeof(int), 1, fp);
    220220       
    221             fros = (int*)malloc(*total_number_of_stations*numSrc*sizeof(int));
    222             lros = (int*)malloc(*total_number_of_stations*numSrc*sizeof(int));
     221            fros = (int*) malloc(*total_number_of_stations*numSrc*sizeof(int));
     222            lros = (int*) malloc(*total_number_of_stations*numSrc*sizeof(int));
    223223     
    224             mytgs0 = (struct tgsrwg*)malloc(*total_number_of_stations*sizeof(struct tgsrwg));
    225             mytgs = (struct tgsrwg*)malloc(*total_number_of_stations*sizeof(struct tgsrwg));
     224            mytgs0 = (struct tgsrwg*) malloc(*total_number_of_stations*sizeof(struct tgsrwg));
     225            mytgs = (struct tgsrwg*) malloc(*total_number_of_stations*sizeof(struct tgsrwg));
    226226
    227227            fread(mytgs0, *total_number_of_stations*sizeof(struct tgsrwg), 1, fp);
    … …  
    229229        else
    230230        {
    231             /* check that the mux files are compatible */     
     231            // Check that the mux files are compatible
    232232            fread(&numsta, sizeof(int), 1, fp);
    233233            if(numsta != *total_number_of_stations)
    … …  
    246246                if (mytgs[j].dt != mytgs0[j].dt)
    247247                {
    248                     fprintf(stderr,"%s has different sampling rate to %s\n",
     248                    fprintf(stderr, "%s has different sampling rate to %s\n",
    249249                    muxFileName,
    250250                    muxFileNameArray[0]);
    … …  
    254254                if (mytgs[j].nt != mytgs0[j].nt)
    255255                {
    256                     fprintf(stderr,"%s has different series length to %s\n",
     256                    fprintf(stderr, "%s has different series length to %s\n",
    257257                    muxFileName,
    258258                    muxFileNameArray[0]);
    … …  
    342342                       verbose);
    343343
    344     //FIXME: The params can be removed after the python interface is modified.
    345     params[0] = (double)total_number_of_stations;
    346     params[1] = (double)delta_t;
    347     params[2] = (double)number_of_time_steps;
    348    
    349344    // Apply rule that an empty permutation file means 'take all stations'
    350345    // We could change this later by passing in None instead of the empty
    … …  
    365360        }
    366361    }
     362   
     363    // The params array is used only for passing data back to Python.
     364    params[0] = (double) number_of_selected_stations;
     365    params[1] = (double) delta_t;
     366    params[2] = (double) number_of_time_steps;
     367   
     368   
    367369   
    368370    // Make array(s) to hold demuxed data for stations given in the
    369371    // permutation file
    370     sts_data = (float**)malloc(number_of_selected_stations*sizeof(float*));
     372    sts_data = (float**) malloc(number_of_selected_stations*sizeof(float*));
    371373    if (sts_data == NULL)
    372374    {
    … …  
    380382    {
    381383        // Initialise sts_data to zero
    382         sts_data[i] = (float*)calloc(len_sts_data, sizeof(float));
     384        sts_data[i] = (float*) calloc(len_sts_data, sizeof(float));
    383385        if (sts_data[i] == NULL)
    384386        {
    … …  
    388390    }
    389391
    390     temp_sts_data = (float*)calloc(len_sts_data, sizeof(float));
    391 
    392     muxData = (float*)malloc(numDataMax*sizeof(float));
     392    temp_sts_data = (float*) calloc(len_sts_data, sizeof(float));
     393
     394    muxData = (float*) malloc(numDataMax*sizeof(float));
    393395   
    394396    /* Loop over all sources */
    … …  
    401403           
    402404   
    403    
    404         /* Read in data block from mux2 file */
     405        // Read in data block from mux2 file
    405406        muxFileName = muxFileNameArray[isrc];
    406407        if((fp = fopen(muxFileName, "r")) == NULL)
    … …  
    427428        //     use ista with mux data
    428429        //     use i with the processed data to be returned         
    429        
    430430        for(i = 0; i < number_of_selected_stations; i++)
    431431        {               
    … …  
    461461            N = number_of_time_steps;
    462462           
    463            
    464             /*
    465             printf("station =  %d, source = %d, fros = %f, lros = %f\n",
    466                    ista,
    467                    isrc,
    468                    (float) fros_per_source[ista],
    469                    (float) lros_per_source[ista]);                 
    470             */
    471                    
    472                    
    473                    
    474463            if (isrc == 0) {
    475464                // Assign values for first source
    … …  
    485474                        printf("Adjusting start time for station %d and source %d",
    486475                               ista, isrc);
    487                         printf(" from %f to %f\n", sts_data[i][N+3], (float)fros_per_source[ista]); 
     476                        printf(" from %f to %f\n",
     477                               sts_data[i][N+3],
     478                               (float) fros_per_source[ista]); 
    488479                    }
    489                     sts_data[i][N+3] = (float)fros_per_source[ista];
     480                    sts_data[i][N+3] = (float) fros_per_source[ista];
    490481                }
    491482               
    492                
    493                 if (sts_data[i][N+4] < (float)lros_per_source[ista]) {         
     483                if (sts_data[i][N+4] < (float) lros_per_source[ista]) {         
    494484                    if (verbose) {
    495485                        printf("Adjusting end time for station %d and source %d",
    496486                               ista, isrc);
    497                         printf(" from %f to %f\n", sts_data[i][N+4], (float)lros_per_source[ista]); 
     487                        printf(" from %f to %f\n",
     488                               sts_data[i][N+4],
     489                               (float) lros_per_source[ista]); 
    498490                    }
    499                     sts_data[i][N+4] = (float)lros_per_source[ista];
     491                    sts_data[i][N+4] = (float) lros_per_source[ista];
    500492                }               
    501493            }
    … …  
    593585    {
    594586        PyErr_SetString(PyExc_ValueError,
    595       "ERROR: Memory for muxFileNameArray could not be allocated.");
     587                        "ERROR: Memory for muxFileNameArray could not be allocated.");
    596588        return NULL;
    597589    }
    … …  
    697689    dimensions[1] = num_ts + POFFSET;
    698690   
    699     pydata = (PyArrayObject*)PyArray_FromDims(2, dimensions, PyArray_DOUBLE);
     691    pydata = (PyArrayObject*) PyArray_FromDims(2, dimensions, PyArray_DOUBLE);
    700692    if(pydata == NULL)
    701693    {
    … …  
    719711                {
    720712                    // This gauge has stopped recording but others are
    721             // still recording
     713                    // still recording
    722714                    *(double*)(pydata->data + i*pydata->strides[0]
    723                                     + time*pydata->strides[1]) =
    724                         0.0;
     715                               + time*pydata->strides[1]) =
     716                      0.0;
    725717                }
    726718                else
Note: See TracChangeset for help on using the changeset viewer.