← Back to branch summary

~ubuntu-branches/ubuntu/wily/eso-midas/wily-proposed

« back to all changes in this revision

Viewing changes to stdred/ccdred/src/ccdmosfit.c

  • Committer: Package Import Robot
  • Author(s): Ole Streicher
  • Date: 2014-04-22 14:44:58 UTC
  • Revision ID: package-import@ubuntu.com-20140422144458-okiwi1assxkkiz39
Tags: upstream-13.09pl1.2+dfsg
Import upstream version 13.09pl1.2+dfsg

Show diffs side-by-side

added added

removed removed

Lines of Context:
stdred/ccdred/src/ccdmosfit.c
 
1
/*===========================================================================
 
2
  Copyright (C) 1995-2011 European Southern Observatory (ESO)
 
3
 
 
4
  This program is free software; you can redistribute it and/or 
 
5
  modify it under the terms of the GNU General Public License as 
 
6
  published by the Free Software Foundation; either version 2 of 
 
7
  the License, or (at your option) any later version.
 
8
 
 
9
  This program is distributed in the hope that it will be useful,
 
10
  but WITHOUT ANY WARRANTY; without even the implied warranty of
 
11
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 
12
  GNU General Public License for more details.
 
13
 
 
14
  You should have received a copy of the GNU General Public 
 
15
  License along with this program; if not, write to the Free 
 
16
  Software Foundation, Inc., 675 Massachusetts Ave, Cambridge, 
 
17
  MA 02139, USA.
 
18
 
 
19
  Correspondence concerning ESO-MIDAS should be addressed as follows:
 
20
        Internet e-mail: midas@eso.org
 
21
        Postal address: European Southern Observatory
 
22
                        Data Management Division 
 
23
                        Karl-Schwarzschild-Strasse 2
 
24
                        D 85748 Garching bei Muenchen 
 
25
                        GERMANY
 
26
===========================================================================*/
 
27
 
 
28
/*+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 
29
.IDENTIFICATION: program ccdmosfit
 
30
.KEYWORDS:       bulk data frames, mosaicing, alignment
 
31
.PURPOSE:        Match the individual subraster elements in the input image.
 
32
                 In order to run this program the user should have created 
 
33
                 the output image and the database file with the CREATE/MOSAIC
 
34
                 command.
 
35
                 The actual algorithm (mo_offset) was developend by Mike Regan
 
36
                 of the Univ. of Maryland
 
37
.VERSION
 
38
 
 
39
 110930         last modif
 
40
-----------------------------------------------------------------------*/
 
41
/*
 
42
 * Define _POSIX_SOURCE to indicate
 
43
 * that this is a POSIX program
 
44
 */
 
45
#define _POSIX_SOURCE 1
 
46
 
 
47
/*
 
48
 * definition of the used functions
 
49
 */
 
50
#include <stdio.h>
 
51
#include <string.h>
 
52
#include <stdlib.h>
 
53
#include <math.h>
 
54
#include <ftoc.h>
 
55
#include <midas_def.h>
 
56
#include <computer.h>
 
57
#include <ccd_def.h>
 
58
 
 
59
/*
 
60
 * define some macros and constants
 
61
 */
 
62
#define NCOL       11
 
63
#define MAXIMS     80
 
64
#define MAXIMSONE  MAXIMS+1
 
65
#define MAXLEV     50                  /* maximum number of contour levels */
 
66
#define MAXPIX     512    /* max frame dimension (X,Y) accessed at once */
 
67
#define MAXSIZ     (MAXPIX * MAXPIX)
 
68
 
 
69
extern void MO_TBLRPAR(), MO_SHIFTS(), MO_FSHIFTS(), MO_SUBALIGN();
 
70
extern void sorti(), MO_ZERO(), MO_FITOFF();
 
71
extern int  MO_M2MATCH(), MO_LINKS(), MO_CLINKS(), MO_FLINKS(), MO_OFFSET();
 
72
 
 
73
/*
 
74
 ++++++++++++++++++++++++++++++++++++++++++++++++++
 
75
 * here starts the code of the function
 
76
 ---------------------------------------------------
 
77
*/
 
78
 
 
79
int main()
 
80
 
 
81
{
 
82
int          uni;
 
83
int          imnoi, imnom, imnoc;
 
84
int          tdb;
 
85
int          sizec;
 
86
int          naxis;
 
87
int          npix[3], npixc[3];
 
88
 
 
89
int          allcoldb, allrowdb;
 
90
int          colo[3];
 
91
int          iaux, i;
 
92
int          imsize[2], nimcol, nimrow;
 
93
int          iav;
 
94
int          minpix;
 
95
int          nimages;
 
96
int          npair;
 
97
int          nshift;
 
98
int          nrsub[2];
 
99
int          ncoldb, nrowdb, nsortdb;
 
100
int          nulo;
 
101
int          nranges, ranges[100];
 
102
int          stat;
 
103
int          tr_area[4];
 
104
int          xyref[2];
 
105
int          verbose;
 
106
int          match;
 
107
int          interp;
 
108
int          inull;
 
109
 
 
110
float        data[3];
 
111
float        usrnul;
 
112
float        rnull;
 
113
 
 
114
double       step[3], start[3];
 
115
double       aostep;
 
116
double       dnull;
 
117
 
 
118
char         *pntri, *pntrc;
 
119
char         *cbuff;
 
120
char         cunit[61];
 
121
char         framei[82], framem[82], framec[62];
 
122
char         tabdb[61];
 
123
char         ident[72];
 
124
char         line[81];
 
125
char         med_area[40], im_area[41], matchlist[41];
 
126
char         output[64];
 
127
 
 
128
static char  x_colo[]     = "X_offpix";
 
129
static char  y_colo[]     = "Y_offpix";
 
130
static char  i_colo[]     = "Offset";
 
131
 
 
132
/* 
 
133
 set up MIDAS environment + enable automatic error abort 
 
134
 */
 
135
 
 
136
SCSPRO("mosfit");
 
137
 
 
138
/*
 
139
 Get the table null value
 
140
 */
 
141
TCMNUL(&inull, &rnull, &dnull);
 
142
MO_NULL = 0.0;
 
143
 
 
144
/* 
 
145
 Get input frame list, tables and output frame
 
146
 --------------------------------------------
 
147
 */
 
148
stat  = SCKGETC("CCDIN",1,80,&iav,framei);             /* input mosaic frame */
 
149
stat  = SCIGET(framei,D_R4_FORMAT, F_IO_MODE, F_IMA_TYPE, MO_DIM2, &naxis, 
 
150
                      npix, start, step, ident, cunit, &pntri, &imnoi);
 
151
 
 
152
stat  = SCKGETC("CCDMSK",1,80,&iav,framem);                    /* input mask */
 
153
if (strncmp(framem,"none",4) == 0)                               /* set flag */
 
154
   imnom = -1;
 
155
else
 
156
   stat  = SCFOPN(framem,D_R4_FORMAT,0,F_IMA_TYPE,&imnom);
 
157
 
 
158
/* 
 
159
 Get database, open it and get info
 
160
 Store the table descriptor into the MO variables
 
161
 */
 
162
stat = SCKGETC("TBLDB",1,60,&iav,tabdb);                   /* database table */
 
163
(void) TCTOPN(tabdb, F_IO_MODE, &tdb);
 
164
(void) TCIGET(tdb, &ncoldb, &nrowdb, &nsortdb, &allcoldb, &allrowdb);
 
165
TCCSER(tdb, x_colo, &colo[0]);
 
166
TCCSER(tdb, y_colo, &colo[1]);
 
167
TCCSER(tdb, i_colo, &colo[2]);
 
168
MO_TBLRPAR(tdb, im_area, med_area);                    /* read the desciptor */
 
169
nimages = MO_NXSUB * MO_NYSUB;
 
170
 
 
171
/*
 
172
 Get the output file 
 
173
 */
 
174
stat = SCKGETC("ccdout",1,60,&iav,framec);              /* final output file */
 
175
 
 
176
/* 
 
177
 get the image section 
 
178
 */
 
179
stat = SCKRDI("TR_SEC",1,4,&iav,tr_area,&uni,&nulo);      /* section include */
 
180
 
 
181
/*
 
182
 Column and row number of reference raster
 
183
 */
 
184
stat   = SCKRDI("NRSUB",1,2,&iav,nrsub,&uni,&nulo);   
 
185
MO_NXRSUB = nrsub[0];
 
186
MO_NYRSUB = nrsub[1];
 
187
if (MO_NXRSUB == 0  || MO_NXRSUB < 1 || MO_NXRSUB > MO_NXSUB)
 
188
   MO_NXRSUB = (MO_NXSUB + 1) / 2;
 
189
if (MO_NYRSUB == 0  || MO_NYRSUB < 1 || MO_NYRSUB > MO_NYSUB)
 
190
   MO_NYRSUB = (MO_NYSUB + 1) / 2;
 
191
 
 
192
/*
 
193
 Get the x and y offset of the reference subraster
 
194
 */
 
195
stat   = SCKRDI("XYREF",1,2,&iav,xyref,&uni,&nulo);   
 
196
MO_XREF = xyref[0];
 
197
MO_YREF = xyref[1];
 
198
 
 
199
/* 
 
200
 get the matching subrasters
 
201
 */
 
202
stat = SCKGETC("EXCLUDE",1,40,&iav,matchlist);      /* subraster match list */
 
203
CGN_UPCOPY(matchlist,matchlist,40);                         /* upper case -> */
 
204
if (strncmp(matchlist,"NONE",4) == 0)      
 
205
   nranges = 0;
 
206
else
 
207
   {
 
208
   cbuff  = (char *) ranges;
 
209
   if ( USRINP( 'i', matchlist, 100, cbuff, &nranges ) != 0 )
 
210
      SCETER(20,"*** FATAL: Error in subraster matching list");
 
211
   if (nranges > 1) sorti(nranges,ranges);
 
212
   }
 
213
 
 
214
/*
 
215
 Get the minimum number of pixels
 
216
 */
 
217
stat   = SCKRDI("MINPIX", 1, 1, &iav, &minpix, &uni, &nulo);   
 
218
 
 
219
/* 
 
220
 Get the size of the output
 
221
 */
 
222
stat   = SCKRDI("OSIZE",1,2,&iav,imsize,&uni,&nulo);    /* size of output */
 
223
nimcol = imsize[0];
 
224
nimrow = imsize[1];
 
225
if ( nimcol != 0 && nimcol > 0 && nimcol >= npix[0])
 
226
   npixc[0] = nimcol;
 
227
else 
 
228
   npixc[0] = npix[0];
 
229
 
 
230
if ( nimrow != 0 && nimrow > 0 && nimrow >= npix[1])
 
231
   npixc[1] = nimrow;
 
232
else 
 
233
   npixc[1] = npix[1];
 
234
 
 
235
/* 
 
236
 Get Null value for the output frame
 
237
 */
 
238
stat = SCKGETC("BLANK",1,20,&iav,output);
 
239
if ((output[0] == '+') && (output[1] == '\0'))
 
240
   iaux = 1;                                    /*  use `last' value as Null */
 
241
else
 
242
  {
 
243
  iav = CGN_CNVT(output,2,1,npixc,&usrnul,&aostep);
 
244
  if (iav < 1)
 
245
     SCETER(19,"*** FATAL: Invalid Null value ... ");
 
246
  MO_BLANK = usrnul;
 
247
  }
 
248
 
 
249
/*
 
250
 get the interpolation
 
251
 */
 
252
stat = SCKGETC("INTERPOL",1,40,&iav,output);           /* section to include */
 
253
CGN_UPSTR(output);                                  /* convert to upper case */
 
254
if (strncmp(output,"NEA",3) == 0)                                /* set flag */
 
255
   interp = MO_BINEAREST;
 
256
else if (strncmp(output,"LIN",3) == 0)
 
257
   interp = MO_BILINEAR;
 
258
else if (strncmp(output,"POLY3",5) == 0) 
 
259
   interp = MO_BIPOLY3;
 
260
else if (strncmp(output,"POLY5",5) == 0) 
 
261
   interp = MO_BIPOLY5;
 
262
else if (strncmp(output,"SPLINE3",7) == 0)
 
263
   interp = MO_BISPLINE3;
 
264
else 
 
265
   interp = 999;
 
266
 
 
267
/*
 
268
 Get the output mode 
 
269
 */
 
270
stat = SCKGETC("VERBOSE",1,3,&iav,output);         /* Get the verbose option */
 
271
CGN_UPSTR(output);                                  /* convert to upper case */
 
272
if (strcmp(output,"YES") == 0)                                   /* set flag */
 
273
   {
 
274
   verbose = 1;
 
275
   strcpy(MO_DEFAULT,"NYFXN");
 
276
   }
 
277
else
 
278
   {
 
279
   verbose = 0;
 
280
   strcpy(MO_DEFAULT,"NYFNN");
 
281
   }
 
282
 
 
283
/* 
 
284
 Allocate memory for output frame and zero it
 
285
 */
 
286
sizec   = npixc[0]*npixc[1];
 
287
stat    = SCFCRE(framec,D_R4_FORMAT,F_O_MODE,F_IMA_TYPE,sizec,&imnoc);
 
288
stat    = SCFMAP(imnoc, F_O_MODE, 1, sizec, &iav, &pntrc);
 
289
if (stat != 0)                                   /* get pointer output frame */
 
290
   SCETER(66,"*** FATAL: Could not allocate virtual memory ...");  
 
291
stat    = SCDWRI(imnoc,"NAXIS",&naxis,1,1,&uni);
 
292
stat    = SCDWRI(imnoc,"NPIX",npix,1,naxis,&uni);
 
293
stat    = SCDWRD(imnoc,"START",start,1,naxis,&uni);
 
294
stat    = SCDWRD(imnoc,"STEP",step,1,naxis,&uni);
 
295
stat    = SCDCOP(imnoi, imnoc, 4, "CUNIT");
 
296
sprintf(ident,"Match of subrasters in input frame %s", framei);
 
297
stat    = SCDWRC(imnoc,"IDENT",1,ident,1,72,&uni);
 
298
 
 
299
MO_ZERO(pntrc, npixc, MO_BLANK);
 
300
(void) SCFPUT(imnoc, 1, npixc[0]*npixc[1], pntrc);
 
301
 
 
302
/*
 
303
 Here we do the real work, that is reading and writing the data
 
304
 First, allocate size of one full strip with noutcols
 
305
 */
 
306
SCTPUT(" ");
 
307
sprintf(line,"Input frame:    %s", framei);
 
308
SCTPUT(line);
 
309
sprintf(line,"Database table: %s", tabdb);
 
310
SCTPUT(line);
 
311
sprintf(line,"Output frame:   %s",framec);
 
312
SCTPUT(line);
 
313
sprintf(line,"Number of subrasters (x,y): %d,%d", MO_NXSUB, MO_NYSUB);
 
314
SCTPUT(line);
 
315
 
 
316
/* 
 
317
  Now we have all parameter input available; time to go to work now
 
318
 */
 
319
if (colo[0] < 0 || colo[1] < 0 )
 
320
   SCETER(4,"*** FATAL: X and/or Y offset column(s) not found");
 
321
 
 
322
MO_FLINKS(tdb, colo, MO_DELTAX, MO_DELTAY, MO_DELTAI, nimages, &nshift);
 
323
if (nshift < nimages) 
 
324
   SCETER(4,"*** FATAL: Fewer shifts than subrasters");
 
325
else
 
326
   MO_FSHIFTS(imnoi, imnoc, MO_DELTAX, MO_DELTAY,
 
327
                            MO_IC1, MO_IC2, MO_IL1, MO_IL2,
 
328
                            MO_OC1, MO_OC2, MO_OL1, MO_OL2);
 
329
/*
 
330
 Get the shift in all overlapping regions
 
331
 */
 
332
match = 1;
 
333
MO_OFFSET(imnoi, imnom, tdb, MO_OC1, MO_OL1, &npair, MO_REOFF, MO_COUNT, 
 
334
          MO_IREF, MO_JREF, minpix, verbose); 
 
335
 
 
336
/*
 
337
 Fit the shifts in the overlapping regions
 
338
 */
 
339
MO_FITOFF(nrowdb, npair, MO_REOFF, MO_COUNT, MO_IREF, MO_JREF, 
 
340
          nranges, ranges, MO_DELTAI);
 
341
 
 
342
/*
 
343
 Do the corrections and align
 
344
 */
 
345
MO_SUBALIGN(imnoi, pntri, imnoc, pntrc, tr_area,
 
346
            MO_IC1, MO_IC2, MO_IL1, MO_IL2,
 
347
            MO_OC1, MO_OC2, MO_OL1, MO_OL2,
 
348
            MO_DELTAX, MO_DELTAY, MO_DELTAI,
 
349
            match, interp, verbose);
 
350
 
 
351
/*
 
352
 Store the data in the data base file and close files. Finally exit
 
353
 */
 
354
(void) SCFPUT(imnoc, 1, npixc[0]*npixc[1], pntrc);
 
355
 
 
356
for (i = 0; i < nrowdb; i++) 
 
357
   {
 
358
   data[0] = MO_DELTAX[i];
 
359
   data[1] = MO_DELTAY[i];
 
360
   data[2] = MO_DELTAI[i];
 
361
   stat = TCRWRR(tdb, i+1, 3, colo, data);
 
362
   }
 
363
 
 
364
(void) TCTCLO(tdb);
 
365
return SCSEPI();
 
366
}
 
367
 
 
368
 
 
369
 
 
370