← Back to branch summary

~ubuntu-branches/debian/jessie/eso-midas/jessie

« back to all changes in this revision

Viewing changes to prim/general/libsrc/zima.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:
prim/general/libsrc/zima.c
 
1
/* @(#)zima.c   19.1 (ES0-DMD) 02/25/03 14:01:10 */
 
2
/*===========================================================================
 
3
  Copyright (C) 1995 European Southern Observatory (ESO)
 
4
 
 
5
  This program is free software; you can redistribute it and/or 
 
6
  modify it under the terms of the GNU General Public License as 
 
7
  published by the Free Software Foundation; either version 2 of 
 
8
  the License, or (at your option) any later version.
 
9
 
 
10
  This program is distributed in the hope that it will be useful,
 
11
  but WITHOUT ANY WARRANTY; without even the implied warranty of
 
12
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 
13
  GNU General Public License for more details.
 
14
 
 
15
  You should have received a copy of the GNU General Public 
 
16
  License along with this program; if not, write to the Free 
 
17
  Software Foundation, Inc., 675 Massachusetss Ave, Cambridge, 
 
18
  MA 02139, USA.
 
19
 
 
20
  Corresponding concerning ESO-MIDAS should be addressed as follows:
 
21
        Internet e-mail: midas@eso.org
 
22
        Postal address: European Southern Observatory
 
23
                        Data Management Division 
 
24
                        Karl-Schwarzschild-Strasse 2
 
25
                        D 85748 Garching bei Muenchen 
 
26
                        GERMANY
 
27
===========================================================================*/
 
28
 
 
29
/*+++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 
30
.COPYRIGHT:  Copyright (c) 1994 European Southern Observatory,
 
31
                                         all rights reserved
 
32
.IDENTIFIER  Czima
 
33
.LANGUAGE    C
 
34
.AUTHOR      K. Banse                   IPG-ESO Garching
 
35
.KEYWORDS    bulk data frame, pixel values
 
36
.PURPOSE     get data values at given real pixel no's. of an image
 
37
.ALGORITHM   a list of real pixel no's. is scanned and corresponding data
 
38
             values of input array are stored in output array
 
39
             bilinear interpolation is used to obtain values for real pixels 
 
40
             not falling exactly in the grid of coordinates for pixel no's. 
 
41
             falling on the upward border lines in x or y, interpolation
 
42
             is done only in x or y, i.e. the border grid points are duplicated
 
43
             to the outside... 
 
44
             calculate also dynamic range on the fly
 
45
.INPUT/OUTPUT
 
46
  call as    Czima( p_in, npix, xpixl, ypixl, ndim, p_out, fmin, fmax )
 
47
  
 
48
  input:
 
49
             float *p_in        input array of data
 
50
             int   *npix        dimension of input data
 
51
             float *xpixl       buffer with real pixel no's. in x
 
52
             float *ypixl       buffer with real pixel no's. in y
 
53
             int   ndim         dimension of XPIXL + YPIXL
 
54
  
 
55
  output:
 
56
             float *p_out       output array of data
 
57
             float *fmin        minimum value of output array
 
58
             float *fmax        max value ...
 
59
 
 
60
.ENVIRONment MIDAS
 
61
 
 
62
.VERSIONS    1.00       940606  from ZIMA.FOR    R.M. van Hees
 
63
.VERSIONS    1.10       941122  fix it , KB
 
64
------------------------------------------------------------*/
 
65
 
 
66
/* Define _POSIX_SOURCE to indicate that this is a POSIX program */
 
67
 
 
68
#define  _POSIX_SOURCE 1
 
69
 
 
70
#include <midas_def.h>
 
71
 
 
72
 
 
73
/*
 
74
 
 
75
*/
 
76
 
 
77
#ifdef __STDC__
 
78
void Czima(float *p_in,int *npix,float *xpixl,float *ypixl,int ndim, 
 
79
           float *p_out,float *fmin,float *fmax)
 
80
 
 
81
#else
 
82
void Czima(p_in,npix,xpixl,ypixl,ndim,p_out,fmin,fmax)
 
83
int   *npix, ndim;
 
84
float *p_in, *p_out, *xpixl, *ypixl, *fmin, *fmax;
 
85
#endif
 
86
 
 
87
{
 
88
register int nn;
 
89
register float cindx;
 
90
 
 
91
int   indxx, indxy, indx1, indx2, indx3, indx4;
 
92
float val, xfrac, yfrac;
 
93
 
 
94
 
 
95
 
 
96
if (*npix <= 1)                         /* only 1-D in Y */
 
97
   {
 
98
   nn = ndim/2;
 
99
   indx1 = ypixl[nn] - 1;
 
100
   if (indx1 < 0) 
 
101
      indx1 = 0;
 
102
   else if (indx1 > (npix[1]-1))
 
103
      indx1 = npix[1] - 1;
 
104
   *fmin = *fmax = p_in[indx1];
 
105
 
 
106
   for (nn=0; nn<ndim; nn++)
 
107
      {
 
108
      cindx = ypixl[nn] - 1.0;          /* move to [0,npix-1] */
 
109
      if (cindx <= 0.0)
 
110
         {
 
111
         indx1 = 0;
 
112
         yfrac = 0.0;
 
113
         }
 
114
      else if (cindx > (float) (npix[1]-1))
 
115
         {
 
116
         indx1 = npix[1] - 1;
 
117
         yfrac = 0.0;
 
118
         }
 
119
      else
 
120
         {
 
121
         indx1 = (int) cindx;
 
122
         yfrac = cindx - indx1;
 
123
         }
 
124
      indx2 = indx1 + 1;
 
125
 
 
126
      if (indx2 > (npix[1]-1))
 
127
         val = p_in[indx1];
 
128
      else
 
129
         val = p_in[indx1] + (yfrac * (p_in[indx2] - p_in[indx1]));
 
130
 
 
131
      if ( *fmin > val )
 
132
         *fmin = val;
 
133
      else if ( *fmax < val )
 
134
         *fmax = val;
 
135
 
 
136
      *p_out++ = val;
 
137
      }
 
138
   }
 
139
 
 
140
else if (npix[1] <= 1)                  /* only 1-D in X */
 
141
   {
 
142
   nn = ndim/2;
 
143
   indx1 = xpixl[nn] - 1;
 
144
   if (indx1 < 0) 
 
145
      indx1 = 0;
 
146
   else if (indx1 > (npix[0]-1))
 
147
      indx1 = npix[0] - 1;
 
148
   *fmin = *fmax = p_in[indx1];
 
149
 
 
150
   for (nn=0; nn<ndim; nn++)
 
151
      {
 
152
      cindx = xpixl[nn] - 1.0;          /* move to [0,npix-1] */
 
153
      if (cindx <= 0.0)
 
154
         {
 
155
         indx1 = 0;
 
156
         xfrac = 0.0;
 
157
         }
 
158
      else if (cindx > (float) (npix[0]-1))
 
159
         {
 
160
         indx1 = *npix - 1;
 
161
         xfrac = 0.0;
 
162
         }
 
163
      else
 
164
         {
 
165
         indx1 = (int) cindx;
 
166
         xfrac = cindx - indx1;
 
167
         }
 
168
      indx2 = indx1 + 1;
 
169
 
 
170
      if (indx2 > (*npix-1))
 
171
         val = p_in[indx1];
 
172
      else
 
173
         val = p_in[indx1] + (xfrac * (p_in[indx2] - p_in[indx1]));
 
174
 
 
175
      if ( *fmin > val )
 
176
         *fmin = val;
 
177
      else if ( *fmax < val )
 
178
         *fmax = val;
 
179
 
 
180
      *p_out++ = val;
 
181
      }
 
182
   }
 
183
 
 
184
else                                            /* 2-D image */
 
185
   {
 
186
   int nsize;
 
187
 
 
188
   nsize = npix[0] * npix[1];
 
189
 
 
190
   nn = ndim/2;
 
191
   indx1 = xpixl[nn] - 1;
 
192
   if (indx1 < 0) 
 
193
      indx1 = 0;
 
194
   else if (indx1 > (npix[0]-1))
 
195
      indx1 = npix[0] - 1;
 
196
   indx2 = ypixl[nn] - 1;
 
197
   if (indx2 < 0) 
 
198
      indx2 = 0;
 
199
   else if (indx2 > (npix[1]-1))
 
200
      indx2 = npix[1] - 1;
 
201
   indx1 += (indx2*npix[0]);
 
202
   *fmin = *fmax = p_in[indx1];
 
203
 
 
204
   for (nn=0; nn<ndim; nn++)
 
205
      {
 
206
      cindx = xpixl[nn] - 1.0;          /* move to [0,npix-1] */
 
207
      if (cindx <= 0.0)
 
208
         {
 
209
         indxx = 0;
 
210
         xfrac = 0.0;
 
211
         }
 
212
      else if (cindx >= (float) (npix[0]-1))
 
213
         {
 
214
         indxx = npix[0] - 1;
 
215
         xfrac = 0.0;
 
216
         }
 
217
      else
 
218
         {
 
219
         indxx = (int) cindx;
 
220
         xfrac = cindx - indxx;
 
221
         }
 
222
      
 
223
      cindx = ypixl[nn] - 1.0;          /* move to [0,npix-1] */
 
224
      if (cindx <= 0.0)
 
225
         {
 
226
         indxy = 0;
 
227
         yfrac = 0.0;
 
228
         }
 
229
      else if (cindx >= (float) (npix[1]-1))
 
230
         {
 
231
         indxy = npix[1] - 1;
 
232
         yfrac = 0.0;
 
233
         }
 
234
      else
 
235
         {
 
236
         indxy = (int) cindx;
 
237
         yfrac = cindx - indxy;
 
238
         }
 
239
 
 
240
      indx1 = indxx + ((*npix) * indxy);
 
241
      indx2 = indx1 + 1;
 
242
      indx3 = indx1 + (*npix);
 
243
      if ((indxx + 1) >= *npix )             /* are we out of grid in x? */
 
244
         {
 
245
         if (indx2 >= nsize)                    /* also out in y? */
 
246
            val = p_in[indx1];
 
247
         else
 
248
            val = p_in[indx1] + ( yfrac * (p_in[indx3] - p_in[indx1]) );
 
249
         }
 
250
      else
 
251
         {
 
252
         if (indx3 >= nsize)                    /* are we out of grid in y? */
 
253
            val = p_in[indx1] + ( xfrac * (p_in[indx2] - p_in[indx1]) );
 
254
         else
 
255
            {
 
256
            indx4 = indx3 + 1;
 
257
            val = p_in[indx1] 
 
258
                  + xfrac * (p_in[indx2] - p_in[indx1])
 
259
                  + yfrac * (p_in[indx3] - p_in[indx1])
 
260
                  + xfrac * yfrac * (p_in[indx1] - p_in[indx2]
 
261
                  - p_in[indx3] + p_in[indx4]);
 
262
            }
 
263
         }
 
264
 
 
265
      if ( *fmin > val )
 
266
         *fmin = val;
 
267
      else if ( *fmax < val )
 
268
         *fmax = val;
 
269
 
 
270
      *p_out++ = val;
 
271
      }
 
272
   }
 
273
}