IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 37463


Ignore:
Timestamp:
Oct 4, 2014, 5:57:37 PM (12 years ago)
Author:
eugene
Message:

add fake images to fakeastro (not quite ready)

Location:
branches/eam_branches/ipp-20140904/Ohana/src/fakeastro
Files:
5 added
13 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/Makefile

    r37449 r37463  
    3030$(SRC)/fakeastro.$(ARCH).o           \
    3131$(SRC)/fakeastro_galaxy.$(ARCH).o    \
     32$(SRC)/StarOps.$(ARCH).o    \
    3233$(SRC)/make_fakestars.$(ARCH).o \
    3334$(SRC)/make_subset.$(ARCH).o \
    … …  
    3738$(SRC)/insert_fakestar.$(ARCH).o \
    3839$(SRC)/gaussian.$(ARCH).o \
     40$(SRC)/fakeastro_images.$(ARCH).o \
     41$(SRC)/load_template_images.$(ARCH).o \
     42$(SRC)/make_fake_images.$(ARCH).o \
     43$(SRC)/get_image_patch.$(ARCH).o \
     44$(SRC)/make_fake_stars.$(ARCH).o \
     45$(SRC)/make_fake_stars_catalog.$(ARCH).o \
     46$(SRC)/save_fake_stars.$(ARCH).o \
     47$(SRC)/match_fake_stars.$(ARCH).o \
    3948$(SRC)/remote_hosts.$(ARCH).o
    4049

  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/include/fakeastro.h

    r37453 r37463  
    1414  int found; // assigned to an object?
    1515} FakeAstro_Stars;
     16
     17typedef struct {
     18  Average  average;
     19  Measure  measure;
     20  Lensing *lensing; // optionally carry out the lensing measurements
     21  int found;
     22} Stars;
    1623
    1724/* used in find_matches, find_matches_refstars */
    … …  
    4148char  *INPUT;
    4249
     50char  *IMAGES_INPUT;
     51char  *CATDIR_INPUT;
     52char  *CATDIR_OUTPUT;
     53
    4354int    PARALLEL;
    4455int    PARALLEL_MANUAL;
    … …  
    5061int    FORCE;
    5162
     63float  RADIUS;
     64
    5265SkyRegion UserPatch;
    5366
    … …  
    5871void          ConfigInit          PROTO((int *argc, char **argv));
    5972void          GetConfig           PROTO((char *config, char *field, char *format, int N, void *ptr));
    60 int           args                PROTO((int argc, char **argv));
    61 int           args_client         PROTO((int argc, char **argv));
     73int           args                PROTO((int *argc, char **argv));
     74int           args_client         PROTO((int *argc, char **argv));
    6275
    6376void          set_db              PROTO((FITS_DB *in));
    … …  
    7083void initialize_client (int argc, char **argv);
    7184
    72 int fakeastro_galaxy (SkyList *skylistInput, HostTable *hosts);
     85int fakeastro_galaxy ();
    7386FakeAstro_Stars *make_fakestars (int Nstars);
    7487int sortStars (FakeAstro_Stars *stars, int Nstars);
    … …  
    96109
    97110int strextend (char *input, char *format,...);
     111
     112int fakeastro_images ();
     113Image *load_template_images (int *nimage);
     114Image *make_fake_images (Image *image, int *nfakeImage);
     115
     116SkyRegion *get_image_patch (Image *image);
     117Stars *make_fake_stars (SkyTable *sky, SkyRegion *patch, Image *image, int *nstars);
     118Stars *make_fake_stars_catalog (Stars *stars, int *nstars, SkyRegion *patch, Catalog *catalog, Image *image);
     119int save_fake_stars (SkyTable *sky, SkyRegion *patch, Image *image, Stars *stars, int Nstars);
     120int match_fake_stars (SkyRegion *region, Stars *stars, unsigned int NstarsIn, Catalog *catalog, Image *image);
     121
     122int InitStar (Stars *star);

  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/ConfigInit.c

    r37453 r37463  
    1919  // if (!ScanConfig (config, "ADDSTAR_RADIUS",         "%lf", 0, &ADDSTAR_RADIUS))   ADDSTAR_RADIUS = 1.0;
    2020
    21   // force CATDIR to be absolute (so parallel mode will work)
    22   GetConfig (config, "CATDIR",                 "%s",  0, CATDIR);
    23   char *tmpcatdir = abspath (CATDIR, DVO_MAX_PATH);
    24   strcpy (CATDIR, tmpcatdir);
    25   free (tmpcatdir);
     21  if (FAKEASTRO_OP == OP_GALAXY) {
     22    // force CATDIR to be absolute (so parallel mode will work)
     23    GetConfig (config, "CATDIR",                 "%s",  0, CATDIR);
     24    char *tmpcatdir = abspath (CATDIR, DVO_MAX_PATH);
     25    strcpy (CATDIR, tmpcatdir);
     26    free (tmpcatdir);
     27   
     28    sprintf (ImageCat, "%s/Images.dat", CATDIR);
     29  }
    2630
    2731  GetConfig (config, "GSCFILE",                "%s",  0, GSCFILE);
    … …  
    3034  GetConfig (config, "PHOTCODE_FILE",          "%s",  0, MasterPhotcodeFile);
    3135
    32   sprintf (ImageCat, "%s/Images.dat", CATDIR);
    33 
    3436  if (!ScanConfig (config, "SKY_DEPTH",         "%d",  0, &SKY_DEPTH)) SKY_DEPTH = 2;
    3537  if (!ScanConfig (config, "SKY_TABLE",         "%s",  0, SKY_TABLE)) SKY_TABLE[0] = 0;
     38
     39  /* set the default search radius */
     40  if (!ScanConfig (config, "ADDSTAR_RADIUS", "%f", 0, &RADIUS)) {
     41    GetConfig (config, "RADIUS", "%f", 0, &RADIUS);
     42  }
     43  if (RADIUS < 0.0001) {
     44    fprintf (stderr, "absurd match radius %f\n", RADIUS);
     45    exit (2);
     46  }
    3647
    3748  if (*CATMODE == 0) strcpy (CATMODE, "RAW");
    3849  if (*CATFORMAT == 0) strcpy (CATFORMAT, "ELIXIR");
    3950
     51  char *CATDIR_CHECK = (FAKEASTRO_OP == OP_GALAXY) ? CATDIR : CATDIR_OUTPUT;
     52
     53  // check for existence of CATDIR
     54  struct stat filestat;
     55  int status = stat (CATDIR_CHECK, &filestat);
     56  if (!FORCE && (status == 0)) {
     57    fprintf (stderr, "directory %s exists, refusing to contaminate\n", CATDIR);
     58    exit (1);
     59  }
     60
    4061  /* update master photcode table if not defined */
    41   sprintf (CatdirPhotcodeFile, "%s/Photcodes.dat", CATDIR);
     62  sprintf (CatdirPhotcodeFile, "%s/Photcodes.dat", CATDIR_CHECK);
    4263  if (!LoadPhotcodes (CatdirPhotcodeFile, MasterPhotcodeFile, TRUE)) {
    4364    fprintf (stderr, "error loading photcode table %s\n", CatdirPhotcodeFile);

  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/args.c

    r37458 r37463  
    33void usage_client (void);
    44
    5 int args (int argc, char **argv) {
     5int args (int *argc, char **argv) {
    66
    77  int N;
    … …  
    1010  FAKEASTRO_OP = OP_NONE;
    1111
    12   if ((N = get_argument (argc, argv, "-images"))) {
    13     remove_argument (N, &argc, argv);
     12  if ((N = get_argument (*argc, argv, "-images"))) {
     13    remove_argument (N, argc, argv);
    1414    FAKEASTRO_OP = OP_IMAGES;
    1515  }
    1616
    17   if ((N = get_argument (argc, argv, "-galaxy"))) {
     17  if ((N = get_argument (*argc, argv, "-galaxy"))) {
    1818    if (FAKEASTRO_OP != OP_NONE) usage();
    19     remove_argument (N, &argc, argv);
     19    remove_argument (N, argc, argv);
    2020    FAKEASTRO_OP = OP_GALAXY;
    2121  }
    … …  
    2828  UserPatch.Dmin = -90;
    2929  UserPatch.Dmax = +90;
    30   if ((N = get_argument (argc, argv, "-region"))) {
    31     remove_argument (N, &argc, argv);
    32     UserPatch.Rmin = atof (argv[N]);
    33     remove_argument (N, &argc, argv);
     30  if ((N = get_argument (*argc, argv, "-region"))) {
     31    remove_argument (N, argc, argv);
     32    UserPatch.Rmin = atof (argv[N]);
     33    remove_argument (N, argc, argv);
    3434    UserPatch.Rmax = atof (argv[N]);
    35     remove_argument (N, &argc, argv);
    36     UserPatch.Dmin = atof (argv[N]);
    37     remove_argument (N, &argc, argv);
     35    remove_argument (N, argc, argv);
     36    UserPatch.Dmin = atof (argv[N]);
     37    remove_argument (N, argc, argv);
    3838    UserPatch.Dmax = atof (argv[N]);
    39     remove_argument (N, &argc, argv);
     39    remove_argument (N, argc, argv);
    4040  }
    41   if ((N = get_argument (argc, argv, "-catalog"))) {
    42     remove_argument (N, &argc, argv);
     41  if ((N = get_argument (*argc, argv, "-catalog"))) {
     42    remove_argument (N, argc, argv);
    4343    UserPatch.Rmin = atof (argv[N]);
    4444    UserPatch.Rmax = UserPatch.Rmin + 0.001;
    45     remove_argument (N, &argc, argv);
     45    remove_argument (N, argc, argv);
    4646    UserPatch.Dmin = atof (argv[N]);
    4747    UserPatch.Dmax = UserPatch.Dmin + 0.001;
    48     remove_argument (N, &argc, argv);
     48    remove_argument (N, argc, argv);
    4949  }
    5050
    … …  
    5252  // based on the state of the SkyTable
    5353  PARALLEL = FALSE;
    54   if ((N = get_argument (argc, argv, "-parallel"))) {
     54  if ((N = get_argument (*argc, argv, "-parallel"))) {
    5555    PARALLEL = TRUE;
    56     remove_argument (N, &argc, argv);
     56    remove_argument (N, argc, argv);
    5757  }
    5858  // this is a test mode : rather than launching the remote jobs and waiting for completion,
    5959  // fakeastro will simply list the remote command and wait for the user to signal completion
    6060  PARALLEL_MANUAL = FALSE;
    61   if ((N = get_argument (argc, argv, "-parallel-manual"))) {
     61  if ((N = get_argument (*argc, argv, "-parallel-manual"))) {
    6262    PARALLEL = TRUE; // -parallel-manual implies -parallel
    6363    PARALLEL_MANUAL = TRUE;
    64     remove_argument (N, &argc, argv);
     64    remove_argument (N, argc, argv);
    6565  }
    6666  // this is a test mode : rather than launching the fakeastro_client jobs remotely, they are
    6767  // run in serial via 'system'
    6868  PARALLEL_SERIAL = FALSE;
    69   if ((N = get_argument (argc, argv, "-parallel-serial"))) {
     69  if ((N = get_argument (*argc, argv, "-parallel-serial"))) {
    7070    if (PARALLEL_MANUAL) {
    7171      fprintf (stderr, "ERROR: cannot mix -parallel-manual and -parallel-serial\n");
    … …  
    7474    PARALLEL = TRUE; // -parallel-serial implies -parallel
    7575    PARALLEL_SERIAL = TRUE;
    76     remove_argument (N, &argc, argv);
     76    remove_argument (N, argc, argv);
    7777  }
    7878
    7979  // MaxDensityUse = FALSE;
    80   // if ((N = get_argument (argc, argv, "-max-density"))) {
    81   //   remove_argument (N, &argc, argv);
     80  // if ((N = get_argument (*argc, argv, "-max-density"))) {
     81  //   remove_argument (N, argc, argv);
    8282  //   MaxDensityValue = atof(argv[N]);
    83   //   remove_argument (N, &argc, argv);
     83  //   remove_argument (N, argc, argv);
    8484  //   MaxDensityUse = TRUE;
    8585  // }
    8686
    8787  FORCE = FALSE;
    88   if ((N = get_argument (argc, argv, "-force"))) {
     88  if ((N = get_argument (*argc, argv, "-force"))) {
    8989    FORCE = TRUE;
    90     remove_argument (N, &argc, argv);
     90    remove_argument (N, argc, argv);
    9191  }
    9292
    9393  VERBOSE = VERBOSE2 = FALSE;
    94   if ((N = get_argument (argc, argv, "-v"))) {
     94  if ((N = get_argument (*argc, argv, "-v"))) {
    9595    VERBOSE = TRUE;
    96     remove_argument (N, &argc, argv);
    97   }
    98   if ((N = get_argument (argc, argv, "-vv"))) {
     96    remove_argument (N, argc, argv);
     97  }
     98  if ((N = get_argument (*argc, argv, "-vv"))) {
    9999    VERBOSE = VERBOSE2 = TRUE;
    100     remove_argument (N, &argc, argv);
     100    remove_argument (N, argc, argv);
    101101  }
    102102
    103103  if (FAKEASTRO_OP == OP_IMAGES) {
    104104    // mandatory arguments to fakeastro -images -input catdir -output -catdir -images images.fits
    105     if ((N = get_argument (argc, argv, "-input"))) {
    106       remove_argument (N, &argc, argv);
     105    if ((N = get_argument (*argc, argv, "-input"))) {
     106      remove_argument (N, argc, argv);
    107107      CATDIR_INPUT = strcreate (argv[N]);
    108       remove_argument (N, &argc, argv);
     108      remove_argument (N, argc, argv);
    109109    } else {
    110110      fprintf (stderr, "missing -input (catdir)\n");
    111111      exit (2);
    112112    }
    113     if ((N = get_argument (argc, argv, "-output"))) {
    114       remove_argument (N, &argc, argv);
     113    if ((N = get_argument (*argc, argv, "-output"))) {
     114      remove_argument (N, argc, argv);
    115115      CATDIR_OUTPUT = strcreate (argv[N]);
    116       remove_argument (N, &argc, argv);
     116      remove_argument (N, argc, argv);
    117117    } else {
    118118      fprintf (stderr, "missing -output (catdir)\n");
    119119      exit (2);
    120120    }
    121     if ((N = get_argument (argc, argv, "-images"))) {
    122       remove_argument (N, &argc, argv);
     121    if ((N = get_argument (*argc, argv, "-input-images"))) {
     122      remove_argument (N, argc, argv);
    123123      IMAGES_INPUT = strcreate (argv[N]);
    124       remove_argument (N, &argc, argv);
     124      remove_argument (N, argc, argv);
    125125    } else {
    126126      fprintf (stderr, "missing -images (images)\n");
    … …  
    129129  }
    130130
    131   if (argc != 1) usage ();
    132131  return TRUE;
    133132}
    134133
    135 int args_client (int argc, char **argv) {
     134int args_client (int *argc, char **argv) {
    136135
    137136  int N;
    … …  
    146145
    147146  HOST_ID = 0;
    148   if ((N = get_argument (argc, argv, "-hostID"))) {
    149     remove_argument (N, &argc, argv);
     147  if ((N = get_argument (*argc, argv, "-hostID"))) {
     148    remove_argument (N, argc, argv);
    150149    HOST_ID = atoi (argv[N]);
    151     remove_argument (N, &argc, argv);
     150    remove_argument (N, argc, argv);
    152151  }
    153152  if (!HOST_ID) usage_client();
    154153
    155154  HOSTDIR = NULL;
    156   if ((N = get_argument (argc, argv, "-hostdir"))) {
    157     remove_argument (N, &argc, argv);
     155  if ((N = get_argument (*argc, argv, "-hostdir"))) {
     156    remove_argument (N, argc, argv);
    158157    HOSTDIR = strcreate (argv[N]);
    159     remove_argument (N, &argc, argv);
     158    remove_argument (N, argc, argv);
    160159  }
    161160  if (!HOSTDIR) usage_client();
    162161
    163162  CPT_FILE = NULL;
    164   if ((N = get_argument (argc, argv, "-cpt"))) {
    165     remove_argument (N, &argc, argv);
     163  if ((N = get_argument (*argc, argv, "-cpt"))) {
     164    remove_argument (N, argc, argv);
    166165    CPT_FILE = strcreate (argv[N]);
    167     remove_argument (N, &argc, argv);
     166    remove_argument (N, argc, argv);
    168167  }
    169168  if (!CPT_FILE) usage_client();
    170169 
    171170  INPUT = NULL;
    172   if ((N = get_argument (argc, argv, "-input"))) {
    173     remove_argument (N, &argc, argv);
     171  if ((N = get_argument (*argc, argv, "-input"))) {
     172    remove_argument (N, argc, argv);
    174173    INPUT = strcreate (argv[N]);
    175     remove_argument (N, &argc, argv);
     174    remove_argument (N, argc, argv);
    176175  }
    177176  if (!INPUT) usage_client();
    178177
    179   if ((N = get_argument (argc, argv, "-images"))) {
    180     remove_argument (N, &argc, argv);
     178  if ((N = get_argument (*argc, argv, "-images"))) {
     179    remove_argument (N, argc, argv);
    181180    FAKEASTRO_OP = OP_IMAGES;
    182181  }
    183   if ((N = get_argument (argc, argv, "-galaxy"))) {
     182  if ((N = get_argument (*argc, argv, "-galaxy"))) {
    184183    if (FAKEASTRO_OP != OP_NONE) usage();
    185     remove_argument (N, &argc, argv);
     184    remove_argument (N, argc, argv);
    186185    FAKEASTRO_OP = OP_GALAXY;
    187186  }
    … …  
    193192  UserPatch.Dmin = -90;
    194193  UserPatch.Dmax = +90;
    195   if ((N = get_argument (argc, argv, "-region"))) {
    196     remove_argument (N, &argc, argv);
    197     UserPatch.Rmin = atof (argv[N]);
    198     remove_argument (N, &argc, argv);
     194  if ((N = get_argument (*argc, argv, "-region"))) {
     195    remove_argument (N, argc, argv);
     196    UserPatch.Rmin = atof (argv[N]);
     197    remove_argument (N, argc, argv);
    199198    UserPatch.Rmax = atof (argv[N]);
    200     remove_argument (N, &argc, argv);
    201     UserPatch.Dmin = atof (argv[N]);
    202     remove_argument (N, &argc, argv);
     199    remove_argument (N, argc, argv);
     200    UserPatch.Dmin = atof (argv[N]);
     201    remove_argument (N, argc, argv);
    203202    UserPatch.Dmax = atof (argv[N]);
    204     remove_argument (N, &argc, argv);
     203    remove_argument (N, argc, argv);
    205204  }
    206   if ((N = get_argument (argc, argv, "-catalog"))) {
    207     remove_argument (N, &argc, argv);
     205  if ((N = get_argument (*argc, argv, "-catalog"))) {
     206    remove_argument (N, argc, argv);
    208207    UserPatch.Rmin = atof (argv[N]);
    209208    UserPatch.Rmax = UserPatch.Rmin + 0.001;
    210     remove_argument (N, &argc, argv);
     209    remove_argument (N, argc, argv);
    211210    UserPatch.Dmin = atof (argv[N]);
    212211    UserPatch.Dmax = UserPatch.Dmin + 0.001;
    213     remove_argument (N, &argc, argv);
     212    remove_argument (N, argc, argv);
    214213  }
    215214
    216215  // MaxDensityUse = FALSE;
    217   // if ((N = get_argument (argc, argv, "-max-density"))) {
    218   //   remove_argument (N, &argc, argv);
     216  // if ((N = get_argument (*argc, argv, "-max-density"))) {
     217  //   remove_argument (N, argc, argv);
    219218  //   MaxDensityValue = atof(argv[N]);
    220   //   remove_argument (N, &argc, argv);
     219  //   remove_argument (N, argc, argv);
    221220  //   MaxDensityUse = TRUE;
    222221  // }
    223222
    224223  VERBOSE = VERBOSE2 = FALSE;
    225   if ((N = get_argument (argc, argv, "-v"))) {
     224  if ((N = get_argument (*argc, argv, "-v"))) {
    226225    VERBOSE = TRUE;
    227     remove_argument (N, &argc, argv);
    228   }
    229   if ((N = get_argument (argc, argv, "-vv"))) {
     226    remove_argument (N, argc, argv);
     227  }
     228  if ((N = get_argument (*argc, argv, "-vv"))) {
    230229    VERBOSE = VERBOSE2 = TRUE;
    231     remove_argument (N, &argc, argv);
    232   }
    233 
    234   if (argc != 1) usage_client ();
     230    remove_argument (N, argc, argv);
     231  }
     232
    235233  return TRUE;
    236234}

  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/fakeastro.c

    r37458 r37463  
    22
    33int main (int argc, char **argv) {
     4
     5  gauss_init (50000);
    46
    57  /* get configuration info, args */
    … …  
    911    case OP_GALAXY:
    1012      /* the object analysis is a separate process iterating over catalogs */
    11       fakeastro_galaxy (skylist, hosts);
     13      fakeastro_galaxy ();
    1214      exit (0);
    1315
    1416    case OP_IMAGES:
    15       // fakeastro_images (skylist);
     17      fakeastro_images ();
    1618      exit (0);
    1719

  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/fakeastro_client.c

    r37449 r37463  
    33int main (int argc, char **argv) {
    44
    5   // need to construct these options with args_loadstarpar...
    6   ConfigInit (&argc, argv);
    7   args_client (argc, argv);
     5  initialize_client (argc, argv);
    86
    97  // client is called with a pointer to the file to be loaded

  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/fakeastro_galaxy.c

    r37458 r37463  
    11# include "fakeastro.h"
    22
    3 int fakeastro_galaxy (SkyList *skylistInput, HostTable *hosts) {
     3int fakeastro_galaxy () {
     4
     5  int n, i;
    46
    57  SkyTable *sky = SkyTableLoadOptimal (CATDIR, SKY_TABLE, GSCFILE, TRUE, SKY_DEPTH, VERBOSE);
    68  SkyTableSetFilenames (sky, CATDIR, "cpt");
    7   SkyList *skylist = SkyListByPatch (sky, -1, &UserPatch);
     9  SkyList *skylistInput = SkyListByPatch (sky, -1, &UserPatch);
    810
    911  // load the list of hosts
    … …  
    1719
    1820    // ensure that the paths are absolute path names
    19     int i;
    2021    for (i = 0; i < hosts->Nhosts; i++) {
    2122      char *tmppath = abspath (hosts->hosts[i].pathname, DVO_MAX_PATH);
    … …  
    2728    init_remote_hosts ();
    2829  }
    29 
    30   int n, i;
    31 
    32   // XXX check for an existing CATDIR and exit?
    33 
    34   /* outline
    35 
    36    * generate N stars (random draw from galaxy mode)
    37    * assign them to catalogs and distribute
    38 
    39    */
    40 
    41   gauss_init (50000);
    4230
    4331  int Nloop = 1;
    … …  
    9381  return TRUE;
    9482}
     83
     84  /* outline
     85
     86   * generate N stars (random draw from galaxy mode)
     87   * assign them to catalogs and distribute
     88
     89   */
     90

  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/fakeastro_images.c

    r37458 r37463  
    11# include "fakeastro.h"
    22
    3 int fakeastro_images (HostTable *hosts) {
     3int fakeastro_images () {
     4
     5  int i, j;
     6
     7  FITS_DB db;
     8
     9  /*** update the image table ***/
     10  snprintf (ImageCat, DVO_MAX_PATH, "%s/Images.dat", CATDIR_OUTPUT);
     11
     12  /* setup image table format and lock */
     13  db.mode    = dvo_catalog_catmode (CATMODE);
     14  db.format  = dvo_catalog_catformat (CATFORMAT);
     15  int status = dvo_image_lock (&db, ImageCat, 3600.0, LCK_XCLD);  // shorter timeout?
     16  if (!status) Shutdown ("ERROR: failure to lock image catalog %s", db.filename);
     17
     18  /* load or create the image table */
     19  if (db.dbstate == LCK_EMPTY) {
     20    if (VERBOSE) fprintf (stderr, "can't find %s, creating a new one\n", ImageCat);
     21    dvo_image_create (&db, GetZeroPoint());
     22  } else {
     23    if (!dvo_image_load (&db, VERBOSE, FALSE)) {
     24      Shutdown ("can't read image catalog %s", db.filename);
     25    }
     26  }
    427
    528  SkyTable *skyTableInput = SkyTableLoadOptimal (CATDIR_INPUT, SKY_TABLE, GSCFILE, TRUE, SKY_DEPTH, VERBOSE);
    629  SkyTableSetFilenames (skyTableInput, CATDIR_INPUT, "cpt");
    7   // SkyList *skyListInput = SkyListByPatch (skyInpt, -1, &UserPatch);
    830
    931  SkyTable *skyTableOutput = SkyTableLoadOptimal (CATDIR_OUTPUT, SKY_TABLE, GSCFILE, TRUE, SKY_DEPTH, VERBOSE);
    1032  SkyTableSetFilenames (skyTableOutput, CATDIR_OUTPUT, "cpt");
    11   // SkyList *skyListOutput = SkyListByPatch (skyInpt, -1, &UserPatch);
    1233
    13   int Nimages;
    14   Images *images = load_template_images (&Nimages);
     34  int Nimage = 0;
     35  int NIMAGE = 1000;
     36  Image *image = NULL;
    1537 
    16   for (i = 0; i < Nimages; i++) {
     38  ALLOCATE (image, Image, NIMAGE);
     39
     40  int Nrefimage;
     41  Image *refimage = load_template_images (&Nrefimage);
     42 
     43  for (i = 0; i < Nrefimage; i++) {
    1744
    1845    // we only want to make fake images for the exposures
    19     if (strcmp(&images[i].coords.ctype[4], "-DIS")) continue;
     46    if (strcmp(&refimage[i].coords.ctype[4], "-DIS")) continue;
    2047
    21     Images *fakeImages = make_fake_images (&images[i], &NfakeImages);
     48    int NfakeImage;
     49    Image *fakeImage = make_fake_images (&refimage[i], &NfakeImage);
    2250   
    23     for (j = 0; j < NfakeImages; j++) {
    24       Stars *fakeStars = make_fake_stars (skyTableInput, skyTableOutput, &fakeImages[i], &NfakeStars);
     51    for (j = 0; j < NfakeImage; j++) {
     52
     53      // we only want to make fake stars for the fake chips
     54      if (strcmp(&fakeImage[j].coords.ctype[4], "-WRP")) continue;
     55
     56      SkyRegion *patch = get_image_patch (&fakeImage[j]);
     57
     58      int NfakeStars;
     59      Stars *fakeStars = make_fake_stars (skyTableInput, patch, &fakeImage[j], &NfakeStars);
    2560     
    26       // save fake stars...
     61      save_fake_stars (skyTableOutput, patch, &fakeImage[j], fakeStars, NfakeStars);
     62
     63      free (fakeStars);
     64      free (patch);
    2765    }
     66
     67    if (Nimage + NfakeImage >= NIMAGE) {
     68      NIMAGE += 1000;
     69      REALLOCATE (image, Image, NIMAGE);
     70    }
     71
     72    for (j = 0; j < NfakeImage; j++) {
     73      memcpy (&image[Nimage], &fakeImage[j], sizeof(Image));
     74    }
     75   
     76    free (fakeImage);
    2877  }
     78
     79  /* add the new image and save */
     80  dvo_image_addrows (&db, image, Nimage);
     81  SetProtect (TRUE);
     82  dvo_image_update (&db, VERBOSE);
     83  SetProtect (FALSE);
     84  dvo_image_unlock (&db); /* unlock? */
     85
     86  exit (0);
    2987}
    3088

  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/initialize.c

    r37453 r37463  
    1212  if (get_argument (argc, argv, "--help")) usage();
    1313
     14  args (&argc, argv);
    1415  ConfigInit (&argc, argv);
    15   args (argc, argv);
    16 
    17   // check for existence of CATDIR
    18   struct stat filestat;
    19   int status = stat (CATDIR, &filestat);
    20   if (!FORCE && (status == 0)) {
    21     fprintf (stderr, "directory %s exist, refusing to contaminate\n", CATDIR);
    22     exit (1);
    23   }
     16  if (argc != 1) usage ();
    2417
    2518}
    … …  
    2720void initialize_client (int argc, char **argv) {
    2821
     22  args_client (&argc, argv);
    2923  ConfigInit (&argc, argv);
    30   args_client (argc, argv);
     24  if (argc != 1) usage_client ();
    3125}
    3226

  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/make_fake_images.c

    r37457 r37463  
    2121    for (iy = 0; iy < 8; iy++) {
    2222     
     23      if ((ix == 0) && (iy == 0)) continue;
     24      if ((ix == 7) && (iy == 0)) continue;
     25      if ((ix == 0) && (iy == 7)) continue;
     26      if ((ix == 7) && (iy == 7)) continue;
     27
    2328      fakeImage[N] = fakeImage[0];
    2429
    … …  
    3035        fakeImage[N].coords.crpix1 = (ix - 3)*4900;
    3136        fakeImage[N].coords.crpix2 = (iy - 3)*4900;
    32         fakeImage[N].coords.pc1_1 = 1.0
    33         fakeImage[N].coords.pc2_2 = 1.0
     37        fakeImage[N].coords.pc1_1 = 1.0;
     38        fakeImage[N].coords.pc2_2 = 1.0;
    3439      } else {
    3540        fakeImage[N].coords.crpix1 = (3 - ix)*4900;
    3641        fakeImage[N].coords.crpix2 = (4 - iy)*4900;
    37         fakeImage[N].coords.pc1_1 = -1.0
    38         fakeImage[N].coords.pc2_2 = -1.0
     42        fakeImage[N].coords.pc1_1 = -1.0;
     43        fakeImage[N].coords.pc2_2 = -1.0;
    3944      }
    4045
    … …  
    5358      fakeImage[N].NY = 4850;
    5459
    55       fakeImage[N].photcode = 10000 + iy*10 + ix;
     60      // XXX need a way to choose this better...
     61      fakeImage[N].photcode = 10100 + iy*10 + ix;
    5662      fakeImage[N].exptime = 45.0;
    5763      N++;

  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/make_fake_stars.c

    r37458 r37463  
    11# include "fakeastro.h"
    22
    3 Stars *make_fake_stars (SkyTable *skyTableInput, SkyTable *skyTableOutput, Image *image, int *nfakeStars) {
     3Stars *make_fake_stars (SkyTable *sky, SkyRegion *patch, Image *image, int *nstars) {
    44
    5   // find the R,D coords of the 4 corners and 4 edge midpoints
    6 
    7   static double Xpt[] = {0, 2425, 4850,    0, 4850,    0, 2425, 4850};
    8   static double Ypt[] = {0,    0,    0, 2425, 2425, 4850, 4850, 4850};
    9 
    10   double Rmin = +480.0;
    11   double Rmax = -360.0;
    12   double Dmin =  +90.0;
    13   double Dmax =  -90.0;
    14 
    15   for (i = 0; i < 8; i++) {
    16     double R, D;
    17     XY_to_RD (&R, &D, Xpt[i], Ypt[i], image->coords);
    18 
    19     Rmin = MIN(R,Rmin);
    20     Rmax = MAX(R,Rmax);
    21     Dmin = MIN(D,Dmin);
    22     Dmax = MAX(D,Dmax);
    23 
    24   }
    25 
    26   SkyList *skyInput  = SkyListByBounds (skyTableInput, -1, Rmin, Rmax, Dmin, Dmax);
     5  SkyList *skylist  = SkyListByPatch (sky, -1, patch);
    276
    287  int Nstars = 0;
     8  Stars *stars = NULL;
     9
     10  Catalog catalog;
    2911
    3012  // load stars from database in these regions
    31   for (i = 0; i < skyInput->Nregions; i++) {
     13  int i;
     14  for (i = 0; i < skylist->Nregions; i++) {
    3215
    3316    // set the parameters which guide catalog open/load/create
    34     catInput.filename  = skyInput[0].filename[i];
    35     catInput.catformat = dvo_catalog_catformat (CATFORMAT);  // set the default catformat from config data
    36     catInput.catmode   = dvo_catalog_catmode (CATMODE);      // set the default catmode from config data
    37     catInput.catflags  = LOAD_AVES | LOAD_SECF;
    38     catInput.Nsecfilt  = GetPhotcodeNsecfilt ();
    39     if (!dvo_catalog_open (&catInput, skylist[0].regions[i], VERBOSE, "r")) {
     17    catalog.filename  = skylist[0].filename[i];
     18    catalog.catformat = dvo_catalog_catformat (CATFORMAT);  // set the default catformat from config data
     19    catalog.catmode   = dvo_catalog_catmode (CATMODE);      // set the default catmode from config data
     20    catalog.catflags  = LOAD_AVES | LOAD_SECF | LOAD_STARPAR;
     21    catalog.Nsecfilt  = GetPhotcodeNsecfilt ();
     22    if (!dvo_catalog_open (&catalog, skylist[0].regions[i], VERBOSE, "r")) {
    4023      continue;
    4124    }
    4225
    4326    // generate fake measurements for this image
    44     stars = make_fake_stars_catalog (catInput, stars, &Nstars);
     27    stars = make_fake_stars_catalog (stars, &Nstars, patch, &catalog, image);
    4528   
     29    dvo_catalog_unlock (&catalog);
     30    dvo_catalog_free (&catalog);
    4631  }
    4732
    48   // will these match, or should I mangle the input cpt names?
    49   SkyList *skyOutput = SkyListByBounds (skyTableOutput, -1, Rmin, Rmax, Dmin, Dmax);
     33  SkyListFree (skylist);
    5034
    51   // load stars from database in these regions
    52   for (i = 0; i < skyOutput->Nregions; i++) {
    53     if (!dvo_catalog_open (&catOutput, skyOutput[0].regions[i], VERBOSE, "w")) {
    54       fprintf (stderr, "ERROR: failure to open/create catalog file %s\n", catalog.filename);
    55       exit (2);
    56     }
    57 
    58     save_fake_stars (stars, Nstars, catOutput);
    59   }
     35  *nstars = Nstars;
     36  return stars;
    6037}

  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/make_fake_stars_catalog.c

    r37460 r37463  
    11# include "fakeastro.h"
     2# define SCALE 0.001
     3
     4// sky[Nsecfilt], crude and hard-wired for now
     5float sky[] = {100.0, 150.0, 200.0, 300.0, 400.0};
    26
    37// things to figure out:
    … …  
    610// * what about QSOs?
    711
    8 Stars *make_fake_stars_catalog (Catalog *catalog, Image *image, int *nfakeStars) {
     12Stars *make_fake_stars_catalog (Stars *stars, int *nstars, SkyRegion *region, Catalog *catalog, Image *image) {
     13
     14  int Nstars = *nstars;
     15  *nstars += catalog->Naverage;
     16
     17  time_t timeRef = ohana_date_to_sec("2011/05/11,00:00:00");
     18
     19  if (!stars) {
     20    ALLOCATE (stars, Stars, *nstars);
     21  } else {
     22    REALLOCATE (stars, Stars, *nstars);
     23  }
     24
     25  int Nsecfilt = GetPhotcodeNsecfilt ();
     26
     27  // use photcode to get zero point
     28  PhotCode *code = GetPhotcodebyCode (image->photcode);
     29  int Nsec = GetPhotcodeNsec (code->equiv);
     30  myAssert (Nsec >= 0, "undefined Nsec?");
    931
    1032  float Mtime = 2.5*log10(image->exptime);
    11   float ZP = image->Mcal + zeropt;
     33
     34  // XXX fix this!!
     35  double plateScale = 0.257;
     36
     37  // XXX put in airmass?
     38  float ZP = SCALE*code->C - image->Mcal + Mtime;
     39  // float ZP = code[0].K*(measure[0].airmass - 1.000) + SCALE*code[0].C - measure[0].Mcal;
    1240
    1341  // generate a set of measurements for each star entry
    1442
    1543  Average *average = catalog->average;
     44  SecFilt *secfilt = catalog->secfilt;
    1645  StarPar *starpar = catalog->starpar;
    1746
     47  int i;
    1848  for (i = 0; i < catalog->Naverage; i++) {
    1949
    20     int nStar = average[i].starparOffset
     50    // true position from src catalog
     51    if (!IN_REGION(average[i].R, average[i].D)) continue;
    2152
    22     InitStar (&stars[i]);
     53    int nStar = average[i].starparOffset;
     54
     55    InitStar (&stars[Nstars]);
    2356
    2457    // which filter?
    25     double Minst = secfilt[i*Nsecfilt + Ns].M  - Mtime - ZP;
     58    double Minst = secfilt[i*Nsecfilt + Nsec].M - ZP;
    2659    double Counts = pow(10.0, -0.4*Minst);
    27     double SkyCts = Something;
     60    double SkyCts = sky[Nsec];
    2861
    2962    double SN = Counts / sqrt(SkyCts + Counts);
    … …  
    3669    // * proper motion
    3770    // * gaussian scatter (~ seeing)
    38     double uR = starpar[nStar].uR;
    39     double uD = starpar[nStar].uD;
     71    double uR = starpar[nStar].uRA / 1000.0; // starpar are (currently) in mas / year
     72    double uD = starpar[nStar].uDEC / 1000.0;
    4073   
    41     double Toffset = average - image.tzero;
     74    // tzero, timeRef are in UNIX seconds, Toffset should be in years
     75    double Toffset = (image->tzero - timeRef) / 365.25 / 86400.0;
    4276
    4377    // uR,uD in linear (arcsec / yr)
    44     double dRoff = uR*Toffset;
    45     double dDoff = uD*Toffset;
     78    double dRpm = uR*Toffset;
     79    double dDpm = uD*Toffset;
    4680
    4781    // uR,uD in linear arcsec
    … …  
    4983    double dDsee = rnd_gauss (0.0, 1.0 / SN);
    5084
     85    double dRoff = (dRpm + dRsee) / 3600.0;
     86    double dDoff = (dDpm + dDsee) / 3600.0;
     87
    5188    double Robs = Rtru + dRoff / cos(Dtru*DEG_RAD);
    5289    double Dobs = Dtru + dDoff;
    5390
    5491    double X, Y;
    55     RD_to_XY (&X, &Y, Robs, Dobs, image->coords);
     92    RD_to_XY (&X, &Y, Robs, Dobs, &image->coords);
     93    if (X < 0) continue;
     94    if (Y < 0) continue;
     95    if (X > image->NX) continue;
     96    if (Y > image->NY) continue;
    5697
    57     stars[i].measure.Xccd       = X;
    58     stars[i].measure.Yccd       = Y;
    59     stars[i].measure.dXccd      = 1.0 / SN / plateScale;
    60     stars[i].measure.dYccd      = 1.0 / SN / plateScale;
     98    stars[Nstars].measure.Xccd       = X;
     99    stars[Nstars].measure.Yccd       = Y;
     100    stars[Nstars].measure.dXccd      = 1.0 / SN / plateScale;
     101    stars[Nstars].measure.dYccd      = 1.0 / SN / plateScale;
    61102
    62     // stars[i].measure.posangle   = ToShortDegrees(ps1data[i].posangle);
    63     // stars[i].measure.pltscale   = ps1data[i].pltscale;
     103    // stars[Nstars].measure.posangle   = ToShortDegrees(ps1data[i].posangle);
     104    // stars[Nstars].measure.pltscale   = ps1data[i].pltscale;
    64105
    65     if ((ps1data[i].M >= 0.0) || isnan(ps1data[i].M)) {
    66       stars[i].measure.M      = NAN;
    67     } else {
    68       stars[i].measure.M      = Minst + ZeroPt;
    69     }
    70     stars[i].measure.dM         = 1.0 / SN;
     106    stars[Nstars].measure.M      = Minst + ZP;
     107    stars[Nstars].measure.dM     = 1.0 / SN;
    71108
    72     // stars[i].measure.dMcal      = ps1data[i].dMcal;
    73     stars[i].measure.Sky        = X?;
    74     stars[i].measure.dSky       = X?;
     109    // stars[Nstars].measure.dMcal      = ps1data[i].dMcal;
     110    stars[Nstars].measure.Sky        = sky[Nsec];
     111    stars[Nstars].measure.dSky       = sqrt(sky[Nsec]);
    75112                       
    76     stars[i].measure.photFlags  = ps1data[i].flags;
    77     stars[i].measure.photFlags2 = ps1data[i].flags2;
     113    stars[Nstars].measure.photFlags  = 0;
     114    stars[Nstars].measure.photFlags2 = 0;
    78115
    79116    // this is may optionally be replaced by the internal sequence (see FilterStars.c)
    80     stars[i].measure.detID      = ps1data[i].detID;
     117    stars[Nstars].measure.detID      = Nstars + 1;
     118   
     119    Nstars ++;
    81120  }
    82121
     122  return stars;
    83123}

  • branches/eam_branches/ipp-20140904/Ohana/src/fakeastro/src/make_fakestars.c

    r37453 r37463  
    3535     */
    3636
    37     double z = rnd_gauss (0.0, Z_GAL);
    38     double r = sqrt(drand48()) * R_GAL;
    39     double Lrad = drand48() * 2 * M_PI;
     37    if (i % 100000 == 0) fprintf (stderr, ".");
     38    double z,r,L,B,R,D,Lrad,Brad;
     39
     40    int inPatch = FALSE;
     41    while (!inPatch) {
     42      z = rnd_gauss (0.0, Z_GAL);
     43      r = sqrt(drand48()) * R_GAL;
     44      Lrad = drand48() * 2 * M_PI;
     45      Brad = atan2(z,r);
     46     
     47      L = Lrad*DEG_RAD;
     48      B = Brad*DEG_RAD;
     49     
     50      ApplyTransform (&R, &D, L, B, transform);
     51      if (R < UserPatch.Rmin) continue;
     52      if (R > UserPatch.Rmax) continue;
     53      if (D < UserPatch.Dmin) continue;
     54      if (D > UserPatch.Dmax) continue;
     55      break;
     56    }
    4057
    4158    // double x = r*cos(L);
    … …  
    4360
    4461    double distance = sqrt (SQ(r) + SQ(z));
    45     double Brad = atan2(z,r);
    4662
    4763    double uL_gal = (A_oort * cos(2.0*Lrad) + B_oort) * cos(Brad) * iFkap;
    … …  
    5975    double Mr = rnd_gauss (11.25, 1.0);
    6076   
    61     double L = Lrad*DEG_RAD;
    62     double B = Brad*DEG_RAD;
    63 
    64     double R, D;
    65     ApplyTransform (&R, &D, L, B, transform);
    66  
    6777    // C1, C2 are from http://arxiv.org/pdf/1306.2945v2.pdf
    6878    double Rrad = R*RAD_DEG;
    … …  
    104114    stars[i].starpar.uDEC = uD;
    105115  }
     116  fprintf (stderr, "\n");
    106117  return stars;
    107118}
Note: See TracChangeset for help on using the changeset viewer.