IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Nov 20, 2006, 2:30:27 AM (20 years ago)
Author:
magnier
Message:

improved tests on source fit

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/psModules/test/objects/tap_pmSourceFitModel.c

    r10095 r10098  
    77#include "pstap.h"
    88
    9 bool fitModels (psRandom *seed, float flux, float radius);
    10 bool fitModelFlux (psRandom *seed, float flux, float radius);
     9bool fitModels (psRandom *seed, float flux, float radius, float sigma);
     10
     11// create a single source (sigma =
     12bool fitModelFlux (psRandom *seed, float flux, float radius, float sigma);
    1113
    1214int main (void)
    … …  
    1517    pmSourceFitModelInit (15, 0.01, 1.0, true);
    1618
    17     // psTraceSetLevel ("psLib.math.psMinimizeLMChi2", 6);
    18 
    1919    plan_tests(240);
    2020
    … …  
    2222    psRandom *seed = psRandomAlloc (PS_RANDOM_TAUS, 0);
    2323
    24     # if (0)
     24    static float radius[] = {3.0, 5.0, 7.0, 10.0, 15.0, 25.0};
     25    static float sigma[] = {1.0, 1.5, 2.0};
     26    static float flux[] = {10000.0, 3000.0, 1000.0, 300.0, 100.0, 30.0, 10.0};
    2527
    26         fitModelFlux (seed, 10000.0);
    27     fitModelFlux (seed, 3000.0);
    28     fitModelFlux (seed, 1000.0);
    29     fitModelFlux (seed, 300.0);
    30     fitModelFlux (seed, 100.0);
    31     fitModelFlux (seed, 30.0);
    32     fitModelFlux (seed, 30.0, 4.0);
    33     fitModelFlux (seed, 30.0, 5.0);
    34     fitModelFlux (seed, 30.0, 7.0);
    35     fitModelFlux (seed, 30.0, 10.0);
    36     fitModelFlux (seed, 30.0, 15.0);
    37     fitModelFlux (seed, 30.0, 20.0);
    38     fitModelFlux (seed, 30.0, 25.0);
    39     fitModelFlux (seed, 30.0, 35.0);
    40     # endif
    41 
    42     fitModels (seed, 100.0,  7.0);
    43     fitModels (seed,  80.0,  7.0);
    44     fitModels (seed,  60.0,  7.0);
    45     fitModels (seed,  40.0,  7.0);
    46     fitModels (seed,  20.0,  7.0);
    47 
    48     fitModels (seed,  30.0,  5.0);
    49     fitModels (seed,  30.0, 10.0);
    50     fitModels (seed,  30.0, 15.0);
     28    for (int i = 0; i < sizeof(sigma)/sizeof(float); i++) {
     29        for (int j = 0; j < sizeof(radius)/sizeof(float); j++) {
     30            for (int k = 0; k < sizeof(flux)/sizeof(float); k++) {
     31                fitModels (seed, flux[k], radius[j], sigma[i]);
     32            }
     33        }
     34    }
    5135
    5236    return exit_status();
    … …  
    5943static psVector *par5 = NULL;
    6044
    61 bool fitModels (psRandom *seed, float flux, float radius)
     45bool fitModels (psRandom *seed, float flux, float radius, float sigma)
    6246{
    6347
    6448    psMemId id = psMemGetId();
    6549
    66     diag("test model fit");
     50    diag("test model fit - flux: %f, radius: %f, sigma: %f", flux, radius, sigma);
    6751
    68     par1 = psVectorAllocEmpty (100, PS_TYPE_F32);
    69     par2 = psVectorAllocEmpty (100, PS_TYPE_F32);
    70     par3 = psVectorAllocEmpty (100, PS_TYPE_F32);
    71     par4 = psVectorAllocEmpty (100, PS_TYPE_F32);
    72     par5 = psVectorAllocEmpty (100, PS_TYPE_F32);
     52    par1 = psVectorAllocEmpty (200, PS_TYPE_F32);
     53    par2 = psVectorAllocEmpty (200, PS_TYPE_F32);
     54    par3 = psVectorAllocEmpty (200, PS_TYPE_F32);
     55    par4 = psVectorAllocEmpty (200, PS_TYPE_F32);
     56    par5 = psVectorAllocEmpty (200, PS_TYPE_F32);
    7357
    74     for (int i = 0; i < 100; i++) {
    75         fitModelFlux (seed, flux, radius);
     58    for (int i = 0; i < 200; i++) {
     59        fitModelFlux (seed, flux, radius, sigma);
    7660    }
     61
     62    float signal = 2*M_PI*sigma*sigma*flux;
     63    float noise = sqrt(signal + 4*M_PI*sigma*sigma*(100 + PS_SQR(5)));
     64    float dMag = noise / signal;
     65    float dPos = sigma * dMag;
     66    diag ("signal: %f, noise: %f, dMag: %f, dPos: %f", signal, noise, dMag, dPos);
    7767
    7868    psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEAN | PS_STAT_SAMPLE_STDEV);
    7969    psVectorStats (stats, par1, NULL, NULL, 0);
    80     fprintf (stderr, "Io: %f +/- %f\n", stats->sampleMean, stats->sampleStdev);
    81     // ok_float_tol (stats->sampleStdev, 0.0, 0.0, "Io ref/fit stdev: %f", stats->sampleStdev);
     70    // fprintf (stderr, "Io: %f +/- %f : %f sig\n", stats->sampleMean, stats->sampleStdev, stats->sampleStdev/dMag);
     71    ok ((stats->sampleStdev/dMag < 2.0), "Io ref/fit stdev: %f : %f sigma", stats->sampleStdev, stats->sampleStdev/dMag);
    8272    psVectorStats (stats, par2, NULL, NULL, 0);
    83     fprintf (stderr, "Xo: %f +/- %f\n", stats->sampleMean, stats->sampleStdev);
     73    // fprintf (stderr, "Xo: %f +/- %f : %f sig\n", stats->sampleMean, stats->sampleStdev, stats->sampleStdev/dPos);
     74    ok ((stats->sampleStdev/dPos < 2.0), "Xo ref/fit stdev: %f : %f sigma", stats->sampleStdev, stats->sampleStdev/dPos);
    8475    // ok_float_tol (stats->sampleStdev, 0.0, 0.0, "Xo ref-fit stdev: %f", stats->sampleStdev);
    8576    psVectorStats (stats, par3, NULL, NULL, 0);
    86     fprintf (stderr, "Yo: %f +/- %f\n", stats->sampleMean, stats->sampleStdev);
     77    // fprintf (stderr, "Yo: %f +/- %f : %f sig\n", stats->sampleMean, stats->sampleStdev, stats->sampleStdev/dPos);
     78    ok ((stats->sampleStdev/dPos < 2.0), "Yo ref/fit stdev: %f : %f sigma", stats->sampleStdev, stats->sampleStdev/dPos);
    8779    // ok_float_tol (stats->sampleStdev, 0.0, 0.0, "Yo ref-fit stdev: %f", stats->sampleStdev);
    8880    psVectorStats (stats, par4, NULL, NULL, 0);
    89     fprintf (stderr, "Sx: %f +/- %f\n", stats->sampleMean, stats->sampleStdev);
     81    // fprintf (stderr, "Sx: %f +/- %f : %f sig\n", stats->sampleMean, stats->sampleStdev, stats->sampleStdev/dMag);
     82    ok ((stats->sampleStdev/dMag < 2.0), "Sx ref/fit stdev: %f : %f sigma", stats->sampleStdev, stats->sampleStdev/dMag);
    9083    // ok_float_tol (stats->sampleStdev, 0.0, 0.0, "Sx ref/fit stdev: %f", stats->sampleStdev);
    9184    psVectorStats (stats, par5, NULL, NULL, 0);
    92     fprintf (stderr, "Sy: %f +/- %f\n", stats->sampleMean, stats->sampleStdev);
     85    // fprintf (stderr, "Sy: %f +/- %f : %f sig\n", stats->sampleMean, stats->sampleStdev, stats->sampleStdev/dMag);
     86    ok ((stats->sampleStdev/dMag < 2.0), "Sy ref/fit stdev: %f : %f sigma", stats->sampleStdev, stats->sampleStdev/dMag);
    9387    // ok_float_tol (stats->sampleStdev, 0.0, 0.0, "Sy ref/fit stdev: %f", stats->sampleStdev);
    9488
    … …  
    9892    psFree (par4);
    9993    psFree (par5);
     94    psFree (stats);
    10095
    10196    ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
    … …  
    10398}
    10499
    105 bool fitModelFlux (psRandom *seed, float flux, float radius)
     100bool fitModelFlux (psRandom *seed, float flux, float radius, float sigma)
    106101{
    107102
    … …  
    122117    source->modelEXT->params->data.F32[2] = 50;
    123118    source->modelEXT->params->data.F32[3] = 50;
    124     source->modelEXT->params->data.F32[4] = 2.0*sqrt(2.0);
    125     source->modelEXT->params->data.F32[5] = 2.0*sqrt(2.0);
     119    source->modelEXT->params->data.F32[4] = 2.0*sqrt(sigma);
     120    source->modelEXT->params->data.F32[5] = 2.0*sqrt(sigma);
    126121    source->modelEXT->params->data.F32[6] = 0;
    127122
    … …  
    176171    psVectorExtend (par5, 100, 1);
    177172
     173    psFree (rnd);
     174    psFree (source);
     175    psFree (guess);
     176
    178177    return true;
    179178}
Note: See TracChangeset for help on using the changeset viewer.