← Back to branch summary

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

« back to all changes in this revision

Viewing changes to stdred/mos/src/mosnorm.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/mos/src/mosnorm.c
 
1
 
 
2
/* mosnorm.c  
 
3
 
 
4
.VERSION
 
5
 090728         last modif
 
6
*/
 
7
 
 
8
 
 
9
#include <tbldef.h>
 
10
#include <midas_def.h>
 
11
#include <ok.h>
 
12
 
 
13
#include <proto_nrutil.h>
 
14
#include <proto_mutil.h>
 
15
 
 
16
#include <stdio.h>
 
17
 
 
18
 
 
19
 
 
20
#define max(A,B) ((A)  > (B) ? (A) : (B))
 
21
#define min(A,B) ((A)  < (B) ? (A) : (B))
 
22
 
 
23
#define MAXSLITS 100
 
24
 
 
25
#define POLYNOM 0
 
26
#define MEDIAN 1
 
27
 
 
28
int Method, Npix[2], Order;
 
29
double Start[2], Step[2];
 
30
 
 
31
#ifdef __STDC__
 
32
void normalize (float[], float[], float[], float[], int[], int[],
 
33
                int[], int);
 
34
void smooth (float[], float[], int, int);
 
35
#else
 
36
void normalize ();
 
37
void smooth ();
 
38
#endif
 
39
 
 
40
 
 
41
 
 
42
int main ()
 
43
 
 
44
{
 
45
  int parm[5];
 
46
  int i, iav;
 
47
  int nax, kun, knul, tid, imno1, imno2;
 
48
  int col, ncol, nslits, acol, arow, nsort, select;
 
49
  int null[3];
 
50
  int slit[MAXSLITS], upper[MAXSLITS], lower[MAXSLITS];
 
51
 
 
52
  int icol[3];
 
53
 
 
54
  char inframe[60], table[60], outframe[60], method[20];
 
55
  float *inimage, *outimage;
 
56
  char unit[64], inst[72], line[80];
 
57
 
 
58
  float moscol[3];
 
59
  float *buff, *buff2;
 
60
 
 
61
/*         get into Midas and get parameters             */
 
62
 
 
63
  SCSPRO ("mosnorm");
 
64
 
 
65
  strcpy(unit,""), strcpy(inst,""); /* Purify */
 
66
 
 
67
  SCKGETC ("IN_A", 1, 60, &iav, inframe);
 
68
  SCKGETC ("IN_B", 1, 60, &iav, table);
 
69
  SCKGETC ("OUT_A", 1, 60, &iav, outframe);
 
70
  SCKGETC ("INPUTC", 1, 20, &iav, method);
 
71
  SCKRDI ("INPUTI", 1, 4, &iav, parm, &kun, &knul);
 
72
 
 
73
  Order = parm[0];
 
74
 
 
75
  SCTPUT ("\n ----------------------- ");
 
76
  sprintf (line, "Input image:         %s ", inframe);
 
77
  SCTPUT (line);
 
78
  sprintf (line, "Input table:         %s ", table);
 
79
  SCTPUT (line);
 
80
  sprintf (line, "Output image:        %s ", outframe);
 
81
  SCTPUT (line);
 
82
  SCTPUT ("\nInput parameters: \n");
 
83
 
 
84
  if (!strncmp (method, "MED", 3) || !strncmp (method, "med", 3))
 
85
    {
 
86
      Method = MEDIAN;
 
87
      sprintf (line, "Smoothing method: Median");
 
88
      SCTPUT (line);
 
89
      sprintf (line, "Window width:        %i\n", Order);
 
90
      SCTPUT (line);
 
91
    }
 
92
  else
 
93
    {
 
94
      Method = POLYNOM;
 
95
      sprintf (line, "Smoothing method: Polynomial");
 
96
      SCTPUT (line);
 
97
      sprintf (line, "Order of fit:        %i\n", Order);
 
98
      SCTPUT (line);
 
99
    }
 
100
 
 
101
/*          map input and output file                              */
 
102
 
 
103
  SCIGET (inframe, D_R4_FORMAT, F_I_MODE, F_IMA_TYPE,
 
104
          2, &nax, Npix, Start, Step, inst, unit, (char **)&inimage, &imno1);
 
105
 
 
106
  SCIPUT (outframe, D_R4_FORMAT, F_O_MODE, F_IMA_TYPE,
 
107
          nax, Npix, Start, Step, inst, unit, (char **)&outimage, &imno2);
 
108
 
 
109
/*           read table MOS                            */
 
110
 
 
111
  TCTOPN (table, F_I_MODE, &tid);
 
112
  TCIGET (tid, &col, &ncol, &nsort, &acol, &arow);
 
113
  TCCSER (tid, ":slit",&icol[0]);
 
114
  TCCSER (tid, ":ystart",&icol[1]);
 
115
  TCCSER (tid, ":yend",&icol[2]);
 
116
 
 
117
  nslits = 0;
 
118
  for (i = 1; i <= ncol; i++)
 
119
    {
 
120
      TCSGET(tid, i, &select);
 
121
      if (select) {
 
122
        TCRRDR (tid, i, 3, icol, moscol, null);
 
123
        slit[nslits] = (int) moscol[0];
 
124
        lower[nslits] =  (int) ((moscol[1]-Start[1])/Step[1]) + 1;
 
125
        upper[nslits] =  (int) ((moscol[2]-Start[1])/Step[1]) + 1;
 
126
        nslits++;
 
127
      }
 
128
    }
 
129
 
 
130
  TCTCLO (tid);
 
131
 
 
132
/*              reserve space for arrays                */
 
133
 
 
134
  buff = (float *) osmmget (Npix[0] * sizeof (float));
 
135
  buff2 = (float *) osmmget (Npix[0] * sizeof (float));
 
136
 
 
137
  normalize (inimage, outimage, buff, buff2, slit, upper, lower, nslits);
 
138
 
 
139
  osmmfree((char *) buff);
 
140
  osmmfree((char *) buff2);
 
141
 
 
142
  /*          epilogue ...                 */
 
143
 
 
144
  SCSEPI ();
 
145
return 0;
 
146
 
 
147
}
 
148
 
 
149
#ifdef __STDC__
 
150
void
 
151
normalize (float inframe[], float outframe[], float rval[],
 
152
           float yval[], int slit[], int upper[], int lower[],
 
153
           int nslits)
 
154
#else
 
155
void
 
156
normalize ( inframe, outframe, rval, yval, slit, upper, lower, nslits)
 
157
  float inframe[],outframe[],rval[],yval[];
 
158
  int   slit[],upper[],lower[],nslits;
 
159
#endif
 
160
 
 
161
{
 
162
 
 
163
  float nrow;
 
164
  int i, nslit, j, k, jj;
 
165
  char line[80];
 
166
 
 
167
  SCTPUT (" slit no. ");
 
168
 
 
169
/* init outframe  to zero                                  */
 
170
 
 
171
  for (j = 0; j < Npix[1]; j++)
 
172
    {
 
173
      for (i = 0; i < Npix[0]; i++)
 
174
        {
 
175
          jj = j * Npix[0] + i;
 
176
          outframe[jj] = 0.0;
 
177
        }
 
178
    }
 
179
 
 
180
 
 
181
/*          loop over all slitlets                         */
 
182
 
 
183
  for (nslit = 0; nslit < nslits ; nslit++) 
 
184
    {
 
185
      sprintf (line, "    %4i", slit[nslit]);
 
186
      SCTPUT (line);
 
187
 
 
188
/*        average all  rows in one slitlet                 */
 
189
 
 
190
      for (j = 0; j < Npix[0]; j++)
 
191
        rval[j] = 0.0;
 
192
      for (j = lower[nslit] - 1; j < upper[nslit]; j++)
 
193
        {
 
194
          for (k = 0; k < Npix[0]; k++)
 
195
            {
 
196
              jj = j * Npix[0] + k;
 
197
              rval[k] = rval[k] + inframe[jj];
 
198
            }
 
199
        }
 
200
 
 
201
      nrow = upper[nslit] - lower[nslit] + 1.0;
 
202
      for (j = 0; j < Npix[0]; j++)
 
203
        rval[j] = rval[j] / nrow;
 
204
 
 
205
      switch (Method)
 
206
        {
 
207
 
 
208
        case MEDIAN:
 
209
 
 
210
          smooth (rval, yval, Npix[0], Order);
 
211
          break;
 
212
 
 
213
        case POLYNOM:
 
214
 
 
215
          fit_poly (rval, yval, 1.0, 1.0, Npix[0], Order);
 
216
          break;
 
217
 
 
218
        }
 
219
 
 
220
/*      now loop through the rows of one slitlet and normalize       */
 
221
 
 
222
      for (j = lower[nslit] -1 ; j < upper[nslit]; j++)
 
223
        {
 
224
          for (k = 0; k < Npix[0]; k++)
 
225
            {
 
226
              jj = j * Npix[0] + k;
 
227
              outframe[jj] = inframe[jj] / yval[k];
 
228
            }
 
229
        }
 
230
    }
 
231
 
 
232
  SCTPUT (" ----------------------- ");
 
233
}
 
234
 
 
235
#ifdef __STDC__
 
236
void
 
237
smooth (float rval[], float yval[], int npix, int order)
 
238
#else
 
239
void
 
240
smooth ( rval, yval, npix, order )
 
241
  float rval[],yval[];
 
242
  int   npix,order;
 
243
#endif
 
244
 
 
245
/*      median filtering of one row: rval in, yval out  */
 
246
 
 
247
{
 
248
 
 
249
  int j, first, last;
 
250
 
 
251
  first = (order - 1) / 2;
 
252
  last = npix - first - 1;
 
253
 
 
254
  for (j = first; j <= last; j++)
 
255
    {
 
256
      yval[j] = pik_median (order, &rval[j - 2]);
 
257
 
 
258
    }
 
259
 
 
260
  for (j = 0; j < first; j++)
 
261
    yval[j] = yval[first];
 
262
  for (j = last + 1; j < npix; j++)
 
263
    yval[j] = yval[last];
 
264
}
 
265