public software.sextractor

[/] [trunk/] [src/] [makeit.c] - Blame information for rev 305

Details | Compare with Previous | View Log

Line No. Rev Author Line
1 2 bertin
/*
2 233 bertin
*                               makeit.c
3 2 bertin
*
4 233 bertin
* Main program.
5 2 bertin
*
6 233 bertin
*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
7 2 bertin
*
8 233 bertin
*       This file part of:      SExtractor
9 2 bertin
*
10 305 bertin
*       Copyright:              (C) 1993-2013 Emmanuel Bertin -- IAP/CNRS/UPMC
11 2 bertin
*
12 233 bertin
*       License:                GNU General Public License
13
*
14
*       SExtractor is free software: you can redistribute it and/or modify
15
*       it under the terms of the GNU General Public License as published by
16
*       the Free Software Foundation, either version 3 of the License, or
17
*       (at your option) any later version.
18
*       SExtractor is distributed in the hope that it will be useful,
19
*       but WITHOUT ANY WARRANTY; without even the implied warranty of
20
*       MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
21
*       GNU General Public License for more details.
22
*       You should have received a copy of the GNU General Public License
23
*       along with SExtractor. If not, see <http://www.gnu.org/licenses/>.
24
*
25 305 bertin
*       Last modified:          13/02/2013
26 233 bertin
*
27
*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
28 2 bertin
 
29
#ifdef HAVE_CONFIG_H
30
#include        "config.h"
31
#endif
32
 
33
#include        <math.h>
34
#include        <stdio.h>
35
#include        <stdlib.h>
36
#include        <string.h>
37
#include        <time.h>
38
 
39
#include        "define.h"
40
#include        "globals.h"
41
#include        "prefs.h"
42
#include        "fits/fitscat.h"
43
#include        "assoc.h"
44
#include        "back.h"
45
#include        "check.h"
46 173 bertin
#include        "fft.h"
47 2 bertin
#include        "field.h"
48
#include        "filter.h"
49
#include        "growth.h"
50
#include        "interpolate.h"
51 173 bertin
#include        "pattern.h"
52 2 bertin
#include        "psf.h"
53 173 bertin
#include        "profit.h"
54 2 bertin
#include        "som.h"
55
#include        "weight.h"
56 8 bertin
#include        "xml.h"
57 2 bertin
 
58 173 bertin
static int              selectext(char *filename);
59
time_t                  thetimet, thetimet2;
60 284 bertin
extern profitstruct     *theprofit,*thedprofit, *thepprofit,*theqprofit;
61 173 bertin
extern char             profname[][32];
62 219 bertin
double                  dtime;
63 11 bertin
 
64 2 bertin
/******************************** makeit *************************************/
65
/*
66
Manage the whole stuff.
67
*/
68
void    makeit()
69
 
70
  {
71
   checkstruct          *check;
72
   picstruct            *dfield, *field,*pffield[MAXFLAG], *wfield,*dwfield;
73
   catstruct            *imacat;
74
   tabstruct            *imatab;
75 173 bertin
   patternstruct        *pattern;
76 2 bertin
   static time_t        thetime1, thetime2;
77 8 bertin
   struct tm            *tm;
78 284 bertin
   unsigned int         modeltype;
79 208 bertin
   int                  nflag[MAXFLAG], nparam2[2],
80 173 bertin
                        i, nok, ntab, next, ntabmax, forcextflag,
81 297 bertin
                        nima0,nima1, nweight0,nweight1, npsf0,npsf1, npat,npat0;
82 2 bertin
 
83 11 bertin
/* Install error logging */
84 16 bertin
  error_installfunc(write_error);
85 11 bertin
 
86 8 bertin
/* Processing start date and time */
87 219 bertin
  dtime = counter_seconds();
88 8 bertin
  thetimet = time(NULL);
89
  tm = localtime(&thetimet);
90
  sprintf(prefs.sdate_start,"%04d-%02d-%02d",
91
        tm->tm_year+1900, tm->tm_mon+1, tm->tm_mday);
92
  sprintf(prefs.stime_start,"%02d:%02d:%02d",
93
        tm->tm_hour, tm->tm_min, tm->tm_sec);
94
 
95
  NFPRINTF(OUTPUT, "");
96
  QPRINTF(OUTPUT, "----- %s %s started on %s at %s with %d thread%s\n\n",
97
                BANNER,
98
                MYVERSION,
99
                prefs.sdate_start,
100
                prefs.stime_start,
101
                prefs.nthreads,
102
                prefs.nthreads>1? "s":"");
103
 
104 2 bertin
/* Initialize globals variables */
105
  initglob();
106
 
107
  NFPRINTF(OUTPUT, "Setting catalog parameters");
108
  readcatparams(prefs.param_name);
109
  useprefs();                   /* update things accor. to prefs parameters */
110
 
111 297 bertin
/* Check if a specific extension should be loaded */
112
  if ((nima0=selectext(prefs.image_name[0])) != RETURN_ERROR)
113
    {
114
    forcextflag = 1;
115
    ntabmax = next = 1;
116
    }
117
  else
118
    forcextflag = 0;
119
 
120
/* Do the same for other data (but do not force single extension mode) */
121
  nima1 = selectext(prefs.image_name[1]);
122
  nweight0 = selectext(prefs.wimage_name[0]);
123
  nweight1 = selectext(prefs.wimage_name[1]);
124
  if (prefs.dpsf_flag)
125
    {
126
    npsf0 = selectext(prefs.psf_name[0]);
127
    npsf1 = selectext(prefs.psf_name[1]);
128
    }
129
  else
130
    npsf0 = selectext(prefs.psf_name[0]);
131
  for (i=0; i<prefs.nfimage_name; i++)
132
    nflag[i] = selectext(prefs.fimage_name[i]);
133
 
134 234 bertin
  if (prefs.psf_flag)
135 2 bertin
    {
136 293 bertin
/*-- Read the first PSF extension to set up stuff such as context parameters */
137 2 bertin
    NFPRINTF(OUTPUT, "Reading PSF information");
138 5 bertin
    if (prefs.dpsf_flag)
139 284 bertin
      {
140 297 bertin
      thedpsf = psf_load(prefs.psf_name[0],nima0<0? 1 :(npsf0<0? 1:npsf0));
141
      thepsf = psf_load(prefs.psf_name[1], nima1<0? 1 :(npsf1<0? 1:npsf1));
142 284 bertin
      }
143
    else
144 297 bertin
      thepsf = psf_load(prefs.psf_name[0], nima0<0? 1 :(npsf0<0? 1:npsf0));
145 2 bertin
 /*-- Need to check things up because of PSF context parameters */
146
    updateparamflags();
147
    useprefs();
148
    }
149
 
150 173 bertin
  if (prefs.prof_flag)
151
    {
152 227 bertin
#ifdef USE_MODEL
153 194 bertin
    fft_init(prefs.nthreads);
154 173 bertin
/* Create profiles at full resolution */
155
    NFPRINTF(OUTPUT, "Preparing profile models");
156 284 bertin
    modeltype = (FLAG(obj2.prof_offset_flux)? MODEL_BACK : MODEL_NONE)
157 244 bertin
        |(FLAG(obj2.prof_dirac_flux)? MODEL_DIRAC : MODEL_NONE)
158
        |(FLAG(obj2.prof_spheroid_flux)?
159
                (FLAG(obj2.prof_spheroid_sersicn)?
160
                        MODEL_SERSIC : MODEL_DEVAUCOULEURS) : MODEL_NONE)
161
        |(FLAG(obj2.prof_disk_flux)? MODEL_EXPONENTIAL : MODEL_NONE)
162
        |(FLAG(obj2.prof_bar_flux)? MODEL_BAR : MODEL_NONE)
163 284 bertin
        |(FLAG(obj2.prof_arms_flux)? MODEL_ARMS : MODEL_NONE);
164
    theprofit = profit_init(thepsf, modeltype);
165 208 bertin
    changecatparamarrays("VECTOR_MODEL", &theprofit->nparam, 1);
166
    changecatparamarrays("VECTOR_MODELERR", &theprofit->nparam, 1);
167
    nparam2[0] = nparam2[1] = theprofit->nparam;
168 209 bertin
    changecatparamarrays("MATRIX_MODELERR", nparam2, 2);
169 284 bertin
    if (prefs.dprof_flag)
170
      thedprofit = profit_init(thedpsf, modeltype);
171 173 bertin
    if (prefs.pattern_flag)
172
      {
173 293 bertin
      npat0 = prefs.prof_disk_patternvectorsize;
174
      if (npat0<prefs.prof_disk_patternmodvectorsize)
175
        npat0 = prefs.prof_disk_patternmodvectorsize;
176
      if (npat0<prefs.prof_disk_patternargvectorsize)
177
        npat0 = prefs.prof_disk_patternargvectorsize;
178 173 bertin
/*---- Do a copy of the original number of pattern components */
179 293 bertin
      prefs.prof_disk_patternncomp = npat0;
180
      pattern = pattern_init(theprofit, prefs.pattern_type, npat0);
181 173 bertin
      if (FLAG(obj2.prof_disk_patternvector))
182
        {
183
        npat = pattern->size[2];
184
        changecatparamarrays("DISK_PATTERN_VECTOR", &npat, 1);
185
        }
186
      if (FLAG(obj2.prof_disk_patternmodvector))
187
        {
188
        npat = pattern->ncomp*pattern->nfreq;
189
        changecatparamarrays("DISK_PATTERNMOD_VECTOR", &npat, 1);
190
        }
191
      if (FLAG(obj2.prof_disk_patternargvector))
192
        {
193
        npat = pattern->ncomp*pattern->nfreq;
194
        changecatparamarrays("DISK_PATTERNARG_VECTOR", &npat, 1);
195
        }
196
      pattern_end(pattern);
197
      }
198
    QPRINTF(OUTPUT, "Fitting model: ");
199
    for (i=0; i<theprofit->nprof; i++)
200
      {
201
      if (i)
202
        QPRINTF(OUTPUT, "+");
203 244 bertin
      QPRINTF(OUTPUT, "%s", theprofit->prof[i]->name);
204 173 bertin
      }
205
    QPRINTF(OUTPUT, "\n");
206 244 bertin
    if (FLAG(obj2.prof_concentration)|FLAG(obj2.prof_concentration))
207 293 bertin
      {
208 244 bertin
      thepprofit = profit_init(thepsf, MODEL_DIRAC);
209
      theqprofit = profit_init(thepsf, MODEL_EXPONENTIAL);
210 293 bertin
      }
211 227 bertin
#else
212
    error(EXIT_FAILURE,
213
                "*Error*: model-fitting is not supported in this build.\n",
214
                        " Please check your configure options");
215
#endif
216 173 bertin
    }
217
 
218 2 bertin
  if (prefs.filter_flag)
219
    {
220
    NFPRINTF(OUTPUT, "Reading detection filter");
221
    getfilter(prefs.filter_name);       /* get the detection filter */
222
    }
223
 
224
  if (FLAG(obj2.sprob))
225
    {
226
    NFPRINTF(OUTPUT, "Initializing Neural Network");
227
    neurinit();
228
    NFPRINTF(OUTPUT, "Reading Neural Network Weights");
229
    getnnw();
230
    }
231
 
232
  if (prefs.somfit_flag)
233
    {
234
     int        margin;
235
 
236
    thesom = som_load(prefs.som_name);
237
    if ((margin=(thesom->inputsize[1]+1)/2) > prefs.cleanmargin)
238
      prefs.cleanmargin = margin;
239
    if (prefs.somfit_vectorsize>thesom->neurdim)
240
      {
241
      prefs.somfit_vectorsize = thesom->neurdim;
242
      sprintf(gstr,"%d", prefs.somfit_vectorsize);
243
      warning("Dimensionality of the SOM-fit vector limited to ", gstr);
244
      }
245
    }
246
 
247
/* Prepare growth-curve buffer */
248
  if (prefs.growth_flag)
249
    initgrowth();
250
 
251 173 bertin
/* Allocate memory for multidimensional catalog parameter arrays */
252
  alloccatparams();
253
  useprefs();
254
 
255 2 bertin
  if (!(imacat = read_cat(prefs.image_name[0])))
256
    error(EXIT_FAILURE, "*Error*: cannot open ", prefs.image_name[0]);
257
  close_cat(imacat);
258
  imatab = imacat->tab;
259 173 bertin
 
260
  if (!forcextflag)
261 2 bertin
    {
262 173 bertin
    ntabmax = imacat->ntab;
263
/*-- Compute the number of valid input extensions */
264
    next = 0;
265
    for (ntab = 0 ; ntab<imacat->ntab; ntab++, imatab = imatab->nexttab)
266
      {
267
/*---- Check for the next valid image extension */
268
      if ((imatab->naxis < 2)
269 2 bertin
        || !strncmp(imatab->xtension, "BINTABLE", 8)
270
        || !strncmp(imatab->xtension, "ASCTABLE", 8))
271 173 bertin
        continue;
272
      next++;
273
      }
274 2 bertin
    }
275 173 bertin
 
276 2 bertin
  thecat.next = next;
277
 
278
/*-- Init the CHECK-images */
279
  if (prefs.check_flag)
280
    {
281
     checkenum  c;
282
 
283
    NFPRINTF(OUTPUT, "Initializing check-image(s)");
284
    for (i=0; i<prefs.ncheck_type; i++)
285 173 bertin
      if ((c=prefs.check_type[i]) != CHECK_NONE)
286
        {
287
        if (prefs.check[c])
288
           error(EXIT_FAILURE,"*Error*: 2 CHECK_IMAGEs cannot have the same ",
289 2 bertin
                        " CHECK_IMAGE_TYPE");
290 173 bertin
        prefs.check[c] = initcheck(prefs.check_name[i], prefs.check_type[i],
291 2 bertin
                        next);
292 173 bertin
        }
293 2 bertin
    }
294
 
295
  NFPRINTF(OUTPUT, "Initializing catalog");
296
  initcat();
297
 
298 8 bertin
/* Initialize XML data */
299 15 bertin
  if (prefs.xml_flag || prefs.cat_type==ASCII_VO)
300 8 bertin
    init_xml(next);
301 2 bertin
 
302
/* Go through all images */
303 297 bertin
  nok = 0;
304 173 bertin
  for (ntab = 0 ; ntab<ntabmax; ntab++, imatab = imatab->nexttab)
305 2 bertin
    {
306
/*--  Check for the next valid image extension */
307 173 bertin
    if (!forcextflag && ((imatab->naxis < 2)
308 2 bertin
        || !strncmp(imatab->xtension, "BINTABLE", 8)
309 173 bertin
        || !strncmp(imatab->xtension, "ASCTABLE", 8)))
310 2 bertin
      continue;
311 297 bertin
 
312 2 bertin
    nok++;
313
 
314
/*-- Initial time measurement*/
315
    time(&thetime1);
316 297 bertin
    thecat.currext = nok;
317 2 bertin
 
318
    dfield = field = wfield = dwfield = NULL;
319
 
320
    if (prefs.dimage_flag)
321
      {
322
/*---- Init the Detection and Measurement-images */
323 173 bertin
      dfield = newfield(prefs.image_name[0], DETECT_FIELD,
324 297 bertin
        nima0<0? ntab:nima0);
325 173 bertin
      field = newfield(prefs.image_name[1], MEASURE_FIELD,
326 297 bertin
        nima1<0? ntab:nima1);
327 2 bertin
      if ((field->width!=dfield->width) || (field->height!=dfield->height))
328
        error(EXIT_FAILURE, "*Error*: Frames have different sizes","");
329
/*---- Prepare interpolation */
330
      if (prefs.dweight_flag && prefs.interp_type[0] == INTERP_ALL)
331
        init_interpolate(dfield, -1, -1);
332
      if (prefs.interp_type[1] == INTERP_ALL)
333
        init_interpolate(field, -1, -1);
334
      }
335
    else
336
      {
337 173 bertin
      field = newfield(prefs.image_name[0], DETECT_FIELD | MEASURE_FIELD,
338 297 bertin
                nima0<0? ntab:nima0);
339 173 bertin
 
340 305 bertin
/*---- Prepare interpolation */
341 2 bertin
      if ((prefs.dweight_flag || prefs.weight_flag)
342
        && prefs.interp_type[0] == INTERP_ALL)
343
      init_interpolate(field, -1, -1);       /* 0.0 or anything else */
344
      }
345
 
346
/*-- Init the WEIGHT-images */
347
    if (prefs.dweight_flag || prefs.weight_flag)
348
      {
349
       weightenum       wtype;
350
       PIXTYPE  interpthresh;
351
 
352
      if (prefs.nweight_type>1)
353
        {
354
/*------ Double-weight-map mode */
355
        if (prefs.weight_type[1] != WEIGHT_NONE)
356
          {
357
/*-------- First: the "measurement" weights */
358
          wfield = newweight(prefs.wimage_name[1],field,prefs.weight_type[1],
359 305 bertin
                (nima1<0 && prefs.image_name[1])?
360
                        ntab : (nweight1<0?1:nweight1));
361 2 bertin
          wtype = prefs.weight_type[1];
362
          interpthresh = prefs.weight_thresh[1];
363
/*-------- Convert the interpolation threshold to variance units */
364
          weight_to_var(wfield, &interpthresh, 1);
365
          wfield->weight_thresh = interpthresh;
366
          if (prefs.interp_type[1] != INTERP_NONE)
367
            init_interpolate(wfield,
368
                prefs.interp_xtimeout[1], prefs.interp_ytimeout[1]);
369
          }
370
/*------ The "detection" weights */
371
        if (prefs.weight_type[0] != WEIGHT_NONE)
372
          {
373
          interpthresh = prefs.weight_thresh[0];
374
          if (prefs.weight_type[0] == WEIGHT_FROMINTERP)
375
            {
376
            dwfield=newweight(prefs.wimage_name[0],wfield,prefs.weight_type[0],
377 297 bertin
                nima0<0? ntab : (nweight0<0? 1 :nweight0));
378 2 bertin
            weight_to_var(wfield, &interpthresh, 1);
379
            }
380
          else
381
            {
382
            dwfield = newweight(prefs.wimage_name[0], dfield?dfield:field,
383 297 bertin
                prefs.weight_type[0], nima0<0? ntab : (nweight0<0?1:nweight0));
384 2 bertin
            weight_to_var(dwfield, &interpthresh, 1);
385
            }
386
          dwfield->weight_thresh = interpthresh;
387
          if (prefs.interp_type[0] != INTERP_NONE)
388
            init_interpolate(dwfield,
389
                prefs.interp_xtimeout[0], prefs.interp_ytimeout[0]);
390
          }
391
        }
392
      else
393
        {
394
/*------ Single-weight-map mode */
395
        wfield = newweight(prefs.wimage_name[0], dfield?dfield:field,
396 297 bertin
                prefs.weight_type[0], nima0<0? ntab : (nweight0<0?1:nweight0));
397 2 bertin
        wtype = prefs.weight_type[0];
398
        interpthresh = prefs.weight_thresh[0];
399
/*------ Convert the interpolation threshold to variance units */
400
        weight_to_var(wfield, &interpthresh, 1);
401
        wfield->weight_thresh = interpthresh;
402
        if (prefs.interp_type[0] != INTERP_NONE)
403
          init_interpolate(wfield,
404
                prefs.interp_xtimeout[0], prefs.interp_ytimeout[0]);
405
        }
406
      }
407
 
408
/*-- Init the FLAG-images */
409
    for (i=0; i<prefs.nimaflag; i++)
410
      {
411 173 bertin
      pffield[i] = newfield(prefs.fimage_name[i], FLAG_FIELD,
412 297 bertin
                nima0<0? ntab : (nflag[i]<0?1:nflag[i]));
413 2 bertin
      if ((pffield[i]->width!=field->width)
414
        || (pffield[i]->height!=field->height))
415
        error(EXIT_FAILURE,
416
        "*Error*: Incompatible FLAG-map size in ", prefs.fimage_name[i]);
417
      }
418
 
419
/*-- Compute background maps for `standard' fields */
420
    QPRINTF(OUTPUT, dfield? "Measurement image:"
421
                        : "Detection+Measurement image: ");
422 243 bertin
    makeback(field, wfield, prefs.wscale_flag[1]);
423 2 bertin
    QPRINTF(OUTPUT, (dfield || (dwfield&&dwfield->flags^INTERP_FIELD))? "(M)   "
424
                "Background: %-10g RMS: %-10g / Threshold: %-10g \n"
425
                : "(M+D) "
426
                "Background: %-10g RMS: %-10g / Threshold: %-10g \n",
427
        field->backmean, field->backsig, (field->flags & DETECT_FIELD)?
428
        field->dthresh: field->thresh);
429
    if (dfield)
430 234 bertin
      {
431 2 bertin
      QPRINTF(OUTPUT, "Detection image: ");
432
      makeback(dfield, dwfield? dwfield
433 243 bertin
                        : (prefs.weight_type[0] == WEIGHT_NONE?NULL:wfield),
434
                prefs.wscale_flag[0]);
435 2 bertin
      QPRINTF(OUTPUT, "(D)   "
436
                "Background: %-10g RMS: %-10g / Threshold: %-10g \n",
437
        dfield->backmean, dfield->backsig, dfield->dthresh);
438
      }
439
    else if (dwfield && dwfield->flags^INTERP_FIELD)
440
      {
441 305 bertin
      makeback(field, dwfield, prefs.wscale_flag[0]);
442 2 bertin
      QPRINTF(OUTPUT, "(D)   "
443
                "Background: %-10g RMS: %-10g / Threshold: %-10g \n",
444
        field->backmean, field->backsig, field->dthresh);
445
      }
446
 
447
/*-- For interpolated weight-maps, copy the background structure */
448
    if (dwfield && dwfield->flags&(INTERP_FIELD|BACKRMS_FIELD))
449
      copyback(dwfield->reffield, dwfield);
450
    if (wfield && wfield->flags&(INTERP_FIELD|BACKRMS_FIELD))
451
      copyback(wfield->reffield, wfield);
452
 
453
/*-- Prepare learn and/or associations */
454
    if (prefs.assoc_flag)
455
      init_assoc(field);                  /* initialize assoc tasks */
456
 
457
/*-- Update the CHECK-images */
458
    if (prefs.check_flag)
459
      for (i=0; i<MAXCHECK; i++)
460
        if ((check=prefs.check[i]))
461
          reinitcheck(field, check);
462
 
463 297 bertin
    if (!forcextflag && nok>1)
464 293 bertin
      {
465
      if (prefs.psf_flag)
466
        {
467
/*------ Read other PSF extensions */
468
        NFPRINTF(OUTPUT, "Reading PSF information");
469
        psf_end(thepsf, thepsfit);
470
        if (prefs.dpsf_flag)
471
          {
472
          psf_end(thedpsf, thedpsfit);
473 297 bertin
          thedpsf = psf_load(prefs.psf_name[0], nok);
474
          thepsf = psf_load(prefs.psf_name[1], nok);
475 293 bertin
          }
476
        else
477 297 bertin
          thepsf = psf_load(prefs.psf_name[0], nok);
478 293 bertin
        }
479
 
480
      if (prefs.prof_flag)
481
        {
482
/*------ Create profiles at full resolution */
483
        profit_end(theprofit);
484
        theprofit = profit_init(thepsf, modeltype);
485
        if (prefs.dprof_flag)
486
          {
487
          profit_end(thedprofit);
488
          thedprofit = profit_init(thedpsf, modeltype);
489
          }
490
        if (prefs.pattern_flag)
491
          {
492
          pattern = pattern_init(theprofit, prefs.pattern_type, npat0);
493
          pattern_end(pattern);
494
          }
495
        if (FLAG(obj2.prof_concentration)|FLAG(obj2.prof_concentration))
496
          {
497
          profit_end(thepprofit);
498
          profit_end(theqprofit);
499
          thepprofit = profit_init(thepsf, MODEL_DIRAC);
500
          theqprofit = profit_init(thepsf, MODEL_EXPONENTIAL);
501
          }
502
        }
503
      }
504
 
505 2 bertin
/*-- Initialize PSF contexts and workspace */
506
    if (prefs.psf_flag)
507
      {
508
      psf_readcontext(thepsf, field);
509 284 bertin
      psf_init();
510 5 bertin
      if (prefs.dpsf_flag)
511
        {
512
        psf_readcontext(thepsf, dfield);
513 284 bertin
        psf_init();
514 5 bertin
        }
515 2 bertin
      }
516
 
517
/*-- Copy field structures to static ones (for catalog info) */
518
    if (dfield)
519
      {
520
      thefield1 = *field;
521
      thefield2 = *dfield;
522
      }
523
    else
524
      thefield1 = thefield2 = *field;
525
 
526
    if (wfield)
527
      {
528
      thewfield1 = *wfield;
529
      thewfield2 = dwfield? *dwfield: *wfield;
530
      }
531
    else if (dwfield)
532
      thewfield2 = *dwfield;
533
 
534
    reinitcat(field);
535
 
536
/*-- Start the extraction pipeline */
537
    NFPRINTF(OUTPUT, "Scanning image");
538
    scanimage(field, dfield, pffield, prefs.nimaflag, wfield, dwfield);
539
 
540
/*-- Finish the current CHECK-image processing */
541
    if (prefs.check_flag)
542
      for (i=0; i<MAXCHECK; i++)
543
        if ((check=prefs.check[i]))
544
          reendcheck(field, check);
545
 
546
/*-- Final time measurements*/
547
    if (time(&thetime2)!=-1)
548
      {
549
      if (!strftime(thecat.ext_date, 12, "%d/%m/%Y", localtime(&thetime2)))
550
        error(EXIT_FAILURE, "*Internal Error*: Date string too long ","");
551
      if (!strftime(thecat.ext_time, 10, "%H:%M:%S", localtime(&thetime2)))
552
        error(EXIT_FAILURE, "*Internal Error*: Time/date string too long ","");
553
      thecat.ext_elapsed = difftime(thetime2, thetime1);
554
      }
555
 
556
    reendcat();
557
 
558 8 bertin
/* Update XML data */
559 297 bertin
    if (prefs.xml_flag || prefs.cat_type==ASCII_VO)
560
      update_xml(&thecat, dfield? dfield:field, field,
561 9 bertin
        dwfield? dwfield:wfield, wfield);
562 8 bertin
 
563
 
564 2 bertin
/*-- Close ASSOC routines */
565
    end_assoc(field);
566
 
567
    for (i=0; i<prefs.nimaflag; i++)
568
      endfield(pffield[i]);
569
    endfield(field);
570
    if (dfield)
571
      endfield(dfield);
572
    if (wfield)
573
      endfield(wfield);
574
    if (dwfield)
575
      endfield(dwfield);
576
 
577 220 bertin
    QPRINTF(OUTPUT, "      Objects: detected %-8d / sextracted %-8d        \n\n",
578 2 bertin
        thecat.ndetect, thecat.ntotal);
579
    }
580
 
581 297 bertin
  if (nok<=0)
582 2 bertin
    error(EXIT_FAILURE, "Not enough valid FITS image extensions in ",
583
        prefs.image_name[0]);
584
  free_cat(&imacat, 1);
585
 
586
  NFPRINTF(OUTPUT, "Closing files");
587
 
588
/* End CHECK-image processing */
589
  if (prefs.check_flag)
590
    for (i=0; i<MAXCHECK; i++)
591
      {
592
      if ((check=prefs.check[i]))
593
        endcheck(check);
594
      prefs.check[i] = NULL;
595
      }
596
 
597
  if (prefs.filter_flag)
598
    endfilter();
599
 
600
  if (prefs.somfit_flag)
601
    som_end(thesom);
602
 
603
  if (prefs.growth_flag)
604
    endgrowth();
605
 
606 227 bertin
#ifdef USE_MODEL
607 173 bertin
  if (prefs.prof_flag)
608
    {
609
    profit_end(theprofit);
610 284 bertin
    if (prefs.dprof_flag)
611
      profit_end(thedprofit);
612 244 bertin
    if (FLAG(obj2.prof_concentration)|FLAG(obj2.prof_concentration))
613
      {
614
      profit_end(thepprofit);
615
      profit_end(theqprofit);
616
      }
617 173 bertin
    fft_end();
618
    }
619 227 bertin
#endif
620 173 bertin
 
621 234 bertin
  if (prefs.psf_flag)
622 284 bertin
    psf_end(thepsf, thepsfit);
623 2 bertin
 
624 5 bertin
  if (prefs.dpsf_flag)
625 284 bertin
    psf_end(thedpsf, thedpsfit);
626 5 bertin
 
627 2 bertin
  if (FLAG(obj2.sprob))
628
    neurclose();
629
 
630 8 bertin
/* Processing end date and time */
631
  thetimet2 = time(NULL);
632
  tm = localtime(&thetimet2);
633
  sprintf(prefs.sdate_end,"%04d-%02d-%02d",
634
        tm->tm_year+1900, tm->tm_mon+1, tm->tm_mday);
635
  sprintf(prefs.stime_end,"%02d:%02d:%02d",
636
        tm->tm_hour, tm->tm_min, tm->tm_sec);
637 219 bertin
  prefs.time_diff = counter_seconds() - dtime;
638 8 bertin
 
639
/* Write XML */
640 17 bertin
  if (prefs.xml_flag)
641
    write_xml(prefs.xml_name);
642 16 bertin
 
643 17 bertin
  endcat((char *)NULL);
644
 
645 15 bertin
  if (prefs.xml_flag || prefs.cat_type==ASCII_VO)
646
    end_xml();
647 8 bertin
 
648 2 bertin
  return;
649
  }
650
 
651
 
652
/******************************** initglob ***********************************/
653
/*
654
Initialize a few global variables
655
*/
656
void    initglob()
657
  {
658
   int  i;
659
 
660
  for (i=0; i<37; i++)
661
    {
662
    ctg[i] = cos(i*PI/18);
663
    stg[i] = sin(i*PI/18);
664
    }
665
 
666
 
667
  return;
668
  }
669
 
670
 
671 173 bertin
/****** selectext ************************************************************
672
PROTO   int selectext(char *filename)
673
PURPOSE Return the user-selected extension number [%d] from the file name.
674
INPUT   Filename character string.
675
OUTPUT  Extension number, or RETURN_ERROR if nos extension specified.
676
NOTES   The bracket and its extension number are removed from the filename if
677
        found.
678
AUTHOR  E. Bertin (IAP)
679
VERSION 08/10/2007
680
 ***/
681
static int      selectext(char *filename)
682
  {
683
   char *bracl,*bracr;
684
   int  next;
685
 
686
  if (filename && (bracl=strrchr(filename, '[')))
687
    {
688
    *bracl = '\0';
689
    if ((bracr=strrchr(bracl+1, ']')))
690
      *bracr = '\0';
691
    next = strtol(bracl+1, NULL, 0);
692
    return next;
693
    }
694
 
695
  return RETURN_ERROR;
696
  }
697
 
698
 
699 16 bertin
/****** write_error ********************************************************
700
PROTO   int     write_error(char *msg1, char *msg2)
701
PURPOSE Manage files in case of a catched error
702
INPUT   a character string,
703
        another character string
704
OUTPUT  RETURN_OK if everything went fine, RETURN_ERROR otherwise.
705
NOTES   -.
706
AUTHOR  E. Bertin (IAP)
707 17 bertin
VERSION 14/07/2006
708 16 bertin
 ***/
709
void    write_error(char *msg1, char *msg2)
710
  {
711
   char                 error[MAXCHAR];
712
 
713
  sprintf(error, "%s%s", msg1,msg2);
714
  if (prefs.xml_flag)
715
    write_xmlerror(prefs.xml_name, error);
716
 
717
/* Also close existing catalog */
718
  endcat(error);
719
 
720
  end_xml();
721
 
722
  return;
723
  }
724
 
725 219 bertin