IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Sep 8, 2008, 6:10:14 PM (18 years ago)
Author:
Paul Price
Message:

Discovered a problem with threaded implementation of pmDarkApply --- the polynomial is altered within pmDarkApplyScan, so it should not be global across threads. Checked and cleaned up the other threaded detrending functions.

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/psModules/src/detrend/pmDark.c

    r19299 r19432  
    222222
    223223    // these calls allocate and save the requested images on the output analysis metadata
    224     psImage *counts = pmReadoutSetAnalysisImage(output->readouts->data[0], PM_READOUT_STACK_ANALYSIS_COUNT, numCols, numRows, PS_TYPE_U16, 0);
     224    psImage *counts = pmReadoutSetAnalysisImage(output->readouts->data[0], PM_READOUT_STACK_ANALYSIS_COUNT,
     225                                                numCols, numRows, PS_TYPE_U16, 0);
    225226    if (!counts) {
    226227        return false;
    227228    }
    228     psImage *sigma = pmReadoutSetAnalysisImage(output->readouts->data[0], PM_READOUT_STACK_ANALYSIS_SIGMA, numCols, numRows, PS_TYPE_F32, NAN);
     229    psImage *sigma = pmReadoutSetAnalysisImage(output->readouts->data[0], PM_READOUT_STACK_ANALYSIS_SIGMA,
     230                                               numCols, numRows, PS_TYPE_F32, NAN);
    229231    if (!sigma) {
    230232        return false;
    231233    }
    232234
    233     psMetadataAddPtr(output->analysis, PS_LIST_TAIL, PM_DARK_ANALYSIS_ORDINATES, PS_DATA_ARRAY | PS_META_REPLACE, "Dark ordinates", ordinates);
    234     psMetadataAddStr(output->analysis, PS_LIST_TAIL, PM_DARK_ANALYSIS_NORM, PS_META_REPLACE, "Dark normalisation", normConcept);
    235 
    236     psMetadataAddPtr(output->analysis, PS_LIST_TAIL, "DARK.VALUES",  PS_DATA_ARRAY  | PS_META_REPLACE, "Dark values", values);
    237     psMetadataAddPtr(output->analysis, PS_LIST_TAIL, "DARK.RO.MASK", PS_DATA_VECTOR | PS_META_REPLACE, "Dark Readout Mask", roMask);
    238     psMetadataAddPtr(output->analysis, PS_LIST_TAIL, "DARK.NORM",    PS_DATA_VECTOR | PS_META_REPLACE, "Dark norm", norm);
    239     psMetadataAddPtr(output->analysis, PS_LIST_TAIL, "DARK.ORDERS",  PS_DATA_VECTOR | PS_META_REPLACE, "Dark orders", orders);
     235    psMetadataAddPtr(output->analysis, PS_LIST_TAIL, PM_DARK_ANALYSIS_ORDINATES,
     236                     PS_DATA_ARRAY | PS_META_REPLACE, "Dark ordinates", ordinates);
     237    psMetadataAddStr(output->analysis, PS_LIST_TAIL, PM_DARK_ANALYSIS_NORM, PS_META_REPLACE,
     238                     "Dark normalisation", normConcept);
     239
     240    psMetadataAddPtr(output->analysis, PS_LIST_TAIL, "DARK.VALUES",  PS_DATA_ARRAY  | PS_META_REPLACE,
     241                     "Dark values", values);
     242    psMetadataAddPtr(output->analysis, PS_LIST_TAIL, "DARK.RO.MASK", PS_DATA_VECTOR | PS_META_REPLACE,
     243                     "Dark Readout Mask", roMask);
     244    psMetadataAddPtr(output->analysis, PS_LIST_TAIL, "DARK.NORM",    PS_DATA_VECTOR | PS_META_REPLACE,
     245                     "Dark norm", norm);
     246    psMetadataAddPtr(output->analysis, PS_LIST_TAIL, "DARK.ORDERS",  PS_DATA_VECTOR | PS_META_REPLACE,
     247                     "Dark orders", orders);
    240248
    241249    for (int i = 0; i < numTerms; i++) {
    … …  
    272280
    273281    // retrieve the required parameter vectors
    274     psArray *values  = psMetadataLookupPtr(&mdok, output->analysis, "DARK.VALUES");  psAssert (values, "values not supplied");
    275     psVector *roMask = psMetadataLookupPtr(&mdok, output->analysis, "DARK.RO.MASK"); psAssert (roMask, "roMask not supplied");
    276     psVector *orders = psMetadataLookupPtr(&mdok, output->analysis, "DARK.ORDERS");  psAssert (orders, "orders not supplied");
     282    psArray *values  = psMetadataLookupPtr(&mdok, output->analysis, "DARK.VALUES");
     283    psAssert(values, "values not supplied");
     284    psVector *roMask = psMetadataLookupPtr(&mdok, output->analysis, "DARK.RO.MASK");
     285    psAssert(roMask, "roMask not supplied");
     286    psVector *orders = psMetadataLookupPtr(&mdok, output->analysis, "DARK.ORDERS");
     287    psAssert(orders, "orders not supplied");
    277288
    278289    psPolynomialMD *poly = psPolynomialMDAlloc(orders); // Polynomial for fitting
    … …  
    307318        int yOut = i - outReadout->row0; // y position on output readout
    308319
    309 # ifdef SHOW_BUSY
     320#ifdef SHOW_BUSY
    310321        if (psTraceGetLevel("psModules.detrend") > 9) {
    311322            printf("Processing row %d\r", i);
    312323            fflush(stdout);
    313324        }
    314 # endif
     325#endif
    315326
    316327        for (int j = minInputCols; j < maxInputCols; j++) {
    … …  
    337348            }
    338349
    339             // XXX test
    340             if (0 && (i == 377) && (j == 80)) {
    341                 FILE *f = fopen ("test.dat", "w");
    342                 for (int r = 0; r < inputs->n; r++) {
    343                     fprintf (f, "%d %d  %d  ", i, j, mask->data.U8[r]);
    344                     psVector *value = values->data[r];
    345                     for (int tmpj = 0; tmpj < value->n; tmpj++) {
    346                         fprintf (f, "%f ", value->data.F32[tmpj]);
    347                     }
    348                     fprintf (f, "%f\n", pixels->data.F32[r]);
    349                 }
    350                 fclose (f);
    351             }
    352 
    353350            if (!psPolynomialMDClipFit(poly, pixels, NULL, mask, maskVal, values, iter, rej)) {
    354351                psErrorClear();         // Nothing we can do about it
    … …  
    367364    psFree(mask);
    368365
    369     fprintf (stderr, "done with combine %x : %f sec\n", (unsigned int) id, psTimerMark (name));
    370366    return true;
    371367}
    … …  
    375371    PS_ASSERT_THREAD_JOB_NON_NULL(job, false);
    376372
    377     pmReadout *readout   = job->args->data[0];
    378     pmCell *dark         = job->args->data[1];
    379     psPolynomialMD *poly = job->args->data[2];
    380     psVector *values     = job->args->data[3];
    381 
    382     psMaskType bad = PS_SCALAR_VALUE(job->args->data[4],U8);
    383     bool doNorm    = PS_SCALAR_VALUE(job->args->data[5],U8);
    384     float norm     = PS_SCALAR_VALUE(job->args->data[6],F32);
    385     int rowStart   = PS_SCALAR_VALUE(job->args->data[7],S32);
    386     int rowStop    = PS_SCALAR_VALUE(job->args->data[8],S32);
    387     bool status = pmDarkApplyScan (readout, dark, poly, values, bad, doNorm, norm, rowStart, rowStop);
    388     return status;
    389 }
    390 
    391 bool pmDarkApplyScan (pmReadout *readout, const pmCell *dark, psPolynomialMD *poly, psVector *values, psMaskType bad, bool doNorm, float norm, int rowStart, int rowStop) {
    392 
     373    pmReadout *readout     = job->args->data[0]; // Readout to correct
     374    pmCell *dark           = job->args->data[1]; // Dark to apply
     375    const psVector *orders = job->args->data[2]; // Polynomial orders for each ordinate
     376    const psVector *values = job->args->data[3]; // Values for each ordinate
     377
     378    psMaskType bad = PS_SCALAR_VALUE(job->args->data[4], U8); // Mask value to give bad pixels
     379    bool doNorm    = PS_SCALAR_VALUE(job->args->data[5], U8); // Normalise values?
     380    float norm     = PS_SCALAR_VALUE(job->args->data[6], F32); // Value by which to normalise
     381    int rowStart   = PS_SCALAR_VALUE(job->args->data[7], S32); // Starting row for scan
     382    int rowStop    = PS_SCALAR_VALUE(job->args->data[8], S32); // Stopping row for scan
     383
     384    return pmDarkApplyScan(readout, dark, orders, values, bad, doNorm, norm, rowStart, rowStop);
     385}
     386
     387bool pmDarkApplyScan(pmReadout *readout, const pmCell *dark, const psVector *orders, const psVector *values,
     388                     psMaskType bad, bool doNorm, float norm, int rowStart, int rowStop)
     389{
    393390    int numCols = readout->image->numCols;
    394391    int numTerms = dark->readouts->n;   // Number of polynomial terms
    395392
    396     // thread here by scan
     393    psPolynomialMD *poly = psPolynomialMDAlloc(orders); // Polynomial to apply
     394
    397395    for (int y = rowStart; y < rowStop; y++) {
    398396        for (int x = 0; x < numCols; x++) {
    … …  
    411409        }
    412410    }
     411    psFree(poly);
     412
    413413    return true;
    414414}
    … …  
    478478    }
    479479
    480     psPolynomialMD *poly = psPolynomialMDAlloc(orders); // Polynomial to apply
    481     psFree(orders);
    482 
    483480    // thread here by scan
    484481
    … …  
    491488
    492489    for (int rowStart = 0; rowStart < readout->image->numRows; rowStart += scanRows) {
    493       int rowStop = PS_MIN (rowStart + scanRows, readout->image->numRows);
    494 
    495       if (threaded) {
    496           // allocate a job, construct the arguments for this job
    497           psThreadJob *job = psThreadJobAlloc ("PSMODULES_DETREND_DARK");
    498           psArrayAdd (job->args, 1, readout);
    499           psArrayAdd (job->args, 1, dark);
    500           psArrayAdd (job->args, 1, poly);
    501           psArrayAdd (job->args, 1, values);
    502           PS_ARRAY_ADD_SCALAR (job->args, bad, PS_TYPE_MASK);
    503           PS_ARRAY_ADD_SCALAR (job->args, doNorm, PS_TYPE_U8);
    504           PS_ARRAY_ADD_SCALAR (job->args, norm, PS_TYPE_F32);
    505           PS_ARRAY_ADD_SCALAR (job->args, rowStart, PS_TYPE_S32);
    506           PS_ARRAY_ADD_SCALAR (job->args, rowStop, PS_TYPE_S32);
    507 
    508           // ppImageDetrendReadout(config, options, view)
    509           if (!psThreadJobAddPending (job)) {
    510               psFree(job);
    511               return false;
    512           }
    513           psFree(job);
    514       } else {
    515           pmDarkApplyScan (readout, dark, poly, values, bad, doNorm, norm, rowStart, rowStop);
    516       }
     490        int rowStop = PS_MIN(rowStart + scanRows, readout->image->numRows);
     491
     492        if (threaded) {
     493            psThreadJob *job = psThreadJobAlloc("PSMODULES_DETREND_DARK");
     494            psArrayAdd(job->args, 1, readout);
     495            psArrayAdd(job->args, 1, dark);
     496            psArrayAdd(job->args, 1, orders);
     497            psArrayAdd(job->args, 1, values);
     498            PS_ARRAY_ADD_SCALAR(job->args, bad, PS_TYPE_MASK);
     499            PS_ARRAY_ADD_SCALAR(job->args, doNorm, PS_TYPE_U8);
     500            PS_ARRAY_ADD_SCALAR(job->args, norm, PS_TYPE_F32);
     501            PS_ARRAY_ADD_SCALAR(job->args, rowStart, PS_TYPE_S32);
     502            PS_ARRAY_ADD_SCALAR(job->args, rowStop, PS_TYPE_S32);
     503
     504            if (!psThreadJobAddPending (job)) {
     505                psFree(job);
     506                psFree(orders);
     507                psFree(values);
     508                return false;
     509            }
     510            psFree(job);
     511        } else if (!pmDarkApplyScan(readout, dark, orders, values, bad, doNorm, norm, rowStart, rowStop)) {
     512            psError(PS_ERR_UNKNOWN, false, "Unable to apply dark.");
     513            psFree(orders);
     514            psFree(values);
     515            return false;
     516        }
    517517    }
    518518
    519519    if (threaded) {
    520520        // wait here for the threaded jobs to finish
    521         if (!psThreadPoolWait(false)) {
    522             psError(PS_ERR_UNKNOWN, false, "Unable to interpolate image.");
     521        if (!psThreadPoolWait(true)) {
     522            psError(PS_ERR_UNKNOWN, false, "Unable to apply dark.");
     523            psFree(orders);
     524            psFree(values);
    523525            return false;
    524526        }
    525         fprintf (stderr, "success for threaded jobs\n");
    526 
    527         // free the done jobs
    528         psThreadJob *job = NULL;
    529         while ((job = psThreadJobGetDone()) != NULL) {
    530             psFree (job);
    531         }
    532     }
    533 
    534     psFree(poly);
     527    }
     528
     529    psFree(orders);
    535530    psFree(values);
    536531
Note: See TracChangeset for help on using the changeset viewer.