source: trunk/MESOSCALE/LMD_MM_MARS/SRC/WPS/geogrid/util/retile.c @ 1243

Last change on this file since 1243 was 11, checked in by aslmd, 14 years ago

spiga@svn-planeto:ajoute le modele meso-echelle martien

File size: 15.5 KB
Line 
1#include <stdio.h>
2#include <sys/types.h>
3#include <sys/stat.h>
4#include <fcntl.h>
5#include <unistd.h>
6#include <math.h>
7#include <string.h>
8
9#define MSG_FLAG 0x80000001
10
11#define N_BYTES 1
12#define IN_TILE_DEGREES_LON 180.0
13#define IN_TILE_DEGREES_LAT 180.0
14#define IN_TILE_PTS_X 1250
15#define IN_TILE_PTS_Y 1250
16#define OUT_TILE_DEGREES_LON 180.0
17#define OUT_TILE_DEGREES_LAT 180.0
18#define OUT_TILE_PTS_X 1250
19#define OUT_TILE_PTS_Y 1250
20#define HALO_WIDTH 3
21#define ZDIM 12
22
23#define CACHESIZE 12
24
25int **** data_cache = NULL;
26char ** fname_cache = NULL;
27int * lru = NULL;
28
29int *** supertile = NULL;
30int supertile_min_x = 999999;
31int supertile_min_y = 999999;
32
33void free_newtile(int *** x)
34{
35  int i, j, k;
36
37  for(i=0; i<IN_TILE_PTS_X; i++)
38  {
39    for(j=0; j<IN_TILE_PTS_Y; j++)
40    {
41      free(x[i][j]);
42    }
43    free(x[i]);
44  }
45  free(x);
46}
47
48int *** read_tile(char * fname)
49{
50  int *** retval;
51  int i, j, k, b;
52  unsigned char buf[IN_TILE_PTS_X*IN_TILE_PTS_Y*ZDIM*N_BYTES];
53  int fd;
54
55  retval = (int ***)malloc(sizeof(int **)*IN_TILE_PTS_X); 
56  for(i=0; i<IN_TILE_PTS_X; i++)
57  {
58    retval[i] = (int **)malloc(sizeof(int *)*IN_TILE_PTS_Y); 
59    for(j=0; j<IN_TILE_PTS_Y; j++)
60    {
61      retval[i][j] = (int *)malloc(sizeof(int)*ZDIM); 
62    }
63  }
64
65  if ((fd = open(fname,O_RDONLY)) == -1) 
66  {
67    fprintf(stderr,"Error opening source file %s\n", fname);
68    return 0;
69  }
70
71  read(fd, (void *)&buf, IN_TILE_PTS_X*IN_TILE_PTS_Y*ZDIM*N_BYTES);
72
73  close(fd);
74
75  /* Put buf into retval */ 
76  for(i=0; i<IN_TILE_PTS_X; i++)
77  {
78    for(j=0; j<IN_TILE_PTS_Y; j++)
79    {
80      for(k=0; k<ZDIM; k++)
81      {
82        retval[i][j][k] = 0;
83        for(b=0; b<N_BYTES; b++)
84        {
85          retval[i][j][k] |= buf[k*N_BYTES*IN_TILE_PTS_X*IN_TILE_PTS_Y+j*N_BYTES*IN_TILE_PTS_X+i*N_BYTES+b] << 8*(N_BYTES-b-1);
86        }
87      }
88    }
89  }
90
91  return retval;
92}
93
94int *** get_tile_from_cache(int i, int j)
95{
96  int ii, jj, kk, k, least, least_idx;
97  int *** retval, *** localptr;
98  char * fname;
99
100  fname = (char *)malloc(256);
101
102  i = (i/IN_TILE_PTS_X)*IN_TILE_PTS_X+1;
103  j = (j/IN_TILE_PTS_Y)*IN_TILE_PTS_Y+1;
104
105  snprintf(fname,256,"%5.5i-%5.5i.%5.5i-%5.5i",i,i+IN_TILE_PTS_X-1,j,j+IN_TILE_PTS_Y-1);
106
107  /* Find out whether tile containing (i,j) is in cache */
108  if (data_cache != NULL) 
109  {
110    for(k=0; k<CACHESIZE; k++)
111    {
112      if (fname_cache[k] != NULL)
113      {
114        if (strncmp(fname_cache[k],fname,256) == 0) 
115        {
116          free(fname);
117          retval = data_cache[k];
118          return retval;
119        }
120      }
121    }
122  }
123
124  /* If not, read from file */
125  localptr = read_tile(fname);
126
127  /* Also store tile in the cache */
128  if (data_cache == NULL)
129  {
130    data_cache = (int ****)malloc(sizeof(int ***)*CACHESIZE);
131    fname_cache = (char **)malloc(sizeof(char *)*CACHESIZE);
132    lru = (int *)malloc(sizeof(int)*CACHESIZE);
133    for(k=0; k<CACHESIZE; k++)
134    {
135      data_cache[k] = NULL;
136      fname_cache[k] = NULL;
137      lru[k] = 0;
138    }
139  }
140
141  least = 0;
142  least_idx = 0;
143  for(k=0; k<CACHESIZE; k++)
144  {
145    lru[k]++; 
146    if (lru[k] > least)
147    {
148      least = lru[k];
149      least_idx = k;
150    }
151  }
152
153  if (data_cache[least_idx] == NULL)
154  {
155    data_cache[least_idx] = localptr;
156    fname_cache[least_idx] = fname; 
157    lru[least_idx] = 0;
158  }
159  else
160  {
161    free_newtile(data_cache[least_idx]);
162    data_cache[least_idx] = localptr;
163    free(fname_cache[least_idx]);
164    fname_cache[least_idx] = fname; 
165    lru[least_idx] = 0;
166  }
167
168  retval = localptr;
169  return retval;
170}
171
172void build_supertile(int i, int j)
173{
174  int ii, jj, kk;
175  int doflip;
176  int *** newtile;
177
178  if (i < 0)
179    i = i + IN_TILE_PTS_X*(int)(360./IN_TILE_DEGREES_LON);
180  if (i >= IN_TILE_PTS_X*(int)(360./IN_TILE_DEGREES_LON))
181    i = i - IN_TILE_PTS_X*(int)(360./IN_TILE_DEGREES_LON);
182  if (j < 0)
183    j = j + IN_TILE_PTS_Y*(int)(180./IN_TILE_DEGREES_LAT);
184  if (j >= IN_TILE_PTS_Y*(int)(180./IN_TILE_DEGREES_LAT))
185    j = j - IN_TILE_PTS_Y*(int)(180./IN_TILE_DEGREES_LAT);
186
187  if (supertile == NULL)
188  {
189    supertile = (int ***)malloc(sizeof(int **)*3*IN_TILE_PTS_X); 
190    for(ii=0; ii<3*IN_TILE_PTS_X; ii++)
191    {
192      supertile[ii] = (int **)malloc(sizeof(int *)*3*IN_TILE_PTS_Y); 
193      for(jj=0; jj<3*IN_TILE_PTS_Y; jj++)
194      {
195        supertile[ii][jj] = (int *)malloc(sizeof(int)*ZDIM); 
196      }
197    }
198  }
199
200  supertile_min_x = (i / IN_TILE_PTS_X)*IN_TILE_PTS_X;
201  supertile_min_y = (j / IN_TILE_PTS_Y)*IN_TILE_PTS_Y;
202
203  /* Get tile containing (i,j) from cache*/
204  /* Get surrounding tiles from cache */ 
205 
206  /* Lower-left */
207  ii = i - IN_TILE_PTS_X;
208  if (ii < 0) ii += (int)(360./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X;
209  doflip = 0;
210  jj = j - IN_TILE_PTS_Y;
211  if (jj < 0)
212  {
213    doflip = 1;
214    jj = j;
215    ii -= (int)(180./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X;
216    if (ii < 0) ii += (int)(360./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X;
217  }
218  newtile = get_tile_from_cache(ii, jj); 
219  if (doflip) 
220  {
221    for(ii=0; ii<IN_TILE_PTS_X; ii++)
222    {
223      for(jj=0; jj<IN_TILE_PTS_Y; jj++)
224      {
225        for(kk=0; kk<ZDIM; kk++)
226          supertile[ii][jj][kk] = newtile[ii][IN_TILE_PTS_Y-1-jj][kk];
227      }
228    }
229  }
230  else
231  {
232    for(ii=0; ii<IN_TILE_PTS_X; ii++)
233    {
234      for(jj=0; jj<IN_TILE_PTS_Y; jj++)
235      {
236        for(kk=0; kk<ZDIM; kk++)
237          supertile[ii][jj][kk] = newtile[ii][jj][kk];
238      }
239    }
240  }
241
242  /* Left */
243  ii = i - IN_TILE_PTS_X;
244  if (ii < 0) ii += (int)(360./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X;
245  jj = j;
246  newtile = get_tile_from_cache(ii, jj); 
247  for(ii=0; ii<IN_TILE_PTS_X; ii++)
248  {
249    for(jj=0; jj<IN_TILE_PTS_Y; jj++)
250    {
251      for(kk=0; kk<ZDIM; kk++)
252        supertile[ii][jj+IN_TILE_PTS_Y][kk] = newtile[ii][jj][kk];
253    }
254  }
255
256  /* Upper-left */
257  ii = i - IN_TILE_PTS_X;
258  if (ii < 0) ii += (int)(360./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X;
259  doflip = 0;
260  jj = j + IN_TILE_PTS_Y;
261  if (jj >= (int)(180./IN_TILE_DEGREES_LAT)*IN_TILE_PTS_Y)
262  {
263    doflip = 1;
264    jj = j;
265    ii -= (int)(180./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X;
266    if (ii < 0) ii += (int)(360./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X;
267  }
268  newtile = get_tile_from_cache(ii, jj); 
269  if (doflip)
270  {
271    for(ii=0; ii<IN_TILE_PTS_X; ii++)
272    {
273      for(jj=0; jj<IN_TILE_PTS_X; jj++)
274      {
275        for(kk=0; kk<ZDIM; kk++)
276          supertile[ii][jj+2*IN_TILE_PTS_Y][kk] = newtile[ii][IN_TILE_PTS_Y-1-jj][kk];
277      }
278    }
279  }
280  else
281  {
282    for(ii=0; ii<IN_TILE_PTS_X; ii++)
283    {
284      for(jj=0; jj<IN_TILE_PTS_Y; jj++)
285      {
286        for(kk=0; kk<ZDIM; kk++)
287          supertile[ii][jj+2*IN_TILE_PTS_Y][kk] = newtile[ii][jj][kk];
288      }
289    }
290  }
291
292  /* Below */
293  ii = i;
294  doflip = 0;
295  jj = j - IN_TILE_PTS_Y;
296  if (jj < 0)
297  {
298    doflip = 1;
299    jj = j;
300    ii -= (int)(180./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X;
301    if (ii < 0) ii += (int)(360./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X;
302  }
303  newtile = get_tile_from_cache(ii, jj); 
304  if (doflip)
305  {
306    for(ii=0; ii<IN_TILE_PTS_X; ii++)
307    {
308      for(jj=0; jj<IN_TILE_PTS_Y; jj++)
309      {
310        for(kk=0; kk<ZDIM; kk++)
311          supertile[ii+IN_TILE_PTS_X][jj][kk] = newtile[ii][IN_TILE_PTS_Y-1-jj][kk];
312      }
313    }
314  }
315  else
316  {
317    for(ii=0; ii<IN_TILE_PTS_X; ii++)
318    {
319      for(jj=0; jj<IN_TILE_PTS_Y; jj++)
320      {
321        for(kk=0; kk<ZDIM; kk++)
322          supertile[ii+IN_TILE_PTS_X][jj][kk] = newtile[ii][jj][kk];
323      }
324    }
325  }
326
327  /* Center */
328  newtile = get_tile_from_cache(i, j); 
329  for(ii=0; ii<IN_TILE_PTS_X; ii++)
330  {
331    for(jj=0; jj<IN_TILE_PTS_Y; jj++)
332    {
333      for(kk=0; kk<ZDIM; kk++)
334        supertile[ii+IN_TILE_PTS_X][jj+IN_TILE_PTS_Y][kk] = newtile[ii][jj][kk];
335    }
336  }
337
338  /* Above */
339  ii = i;
340  doflip = 0;
341  jj = j + IN_TILE_PTS_Y;
342  if (jj >= (int)(180./IN_TILE_DEGREES_LAT)*IN_TILE_PTS_Y)
343  {
344    doflip = 1;
345    jj = j;
346    ii -= (int)(180./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X;
347    if (ii < 0) ii += (int)(360./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X;
348  }
349  newtile = get_tile_from_cache(ii, jj); 
350  if (doflip) 
351  {
352    for(ii=0; ii<IN_TILE_PTS_X; ii++)
353    {
354      for(jj=0; jj<IN_TILE_PTS_Y; jj++)
355      {
356        for(kk=0; kk<ZDIM; kk++)
357          supertile[ii+IN_TILE_PTS_X][jj+2*IN_TILE_PTS_Y][kk] = newtile[ii][IN_TILE_PTS_Y-1-jj][kk];
358      }
359    }
360  }
361  else
362  {
363    for(ii=0; ii<IN_TILE_PTS_X; ii++)
364    {
365      for(jj=0; jj<IN_TILE_PTS_Y; jj++)
366      {
367        for(kk=0; kk<ZDIM; kk++)
368          supertile[ii+IN_TILE_PTS_X][jj+2*IN_TILE_PTS_Y][kk] = newtile[ii][jj][kk];
369      }
370    }
371  }
372
373  /* Lower-right */
374  ii = i + IN_TILE_PTS_X;
375  if (ii >= (int)(360./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X) ii -= (int)(360./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X;
376  doflip = 0;
377  jj = j - IN_TILE_PTS_Y;
378  if (jj < 0)
379  {
380    doflip = 1;
381    jj = j;
382    ii -= (int)(180./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X;
383    if (ii < 0) ii += (int)(360./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X;
384  }
385  newtile = get_tile_from_cache(ii, jj); 
386  if (doflip) 
387  {
388    for(ii=0; ii<IN_TILE_PTS_X; ii++)
389    {
390      for(jj=0; jj<IN_TILE_PTS_Y; jj++)
391      {
392        for(kk=0; kk<ZDIM; kk++)
393          supertile[ii+2*IN_TILE_PTS_X][jj][kk] = newtile[ii][IN_TILE_PTS_Y-1-jj][kk];
394      }
395    }
396  }
397  else
398  {
399    for(ii=0; ii<IN_TILE_PTS_X; ii++)
400    {
401      for(jj=0; jj<IN_TILE_PTS_Y; jj++)
402      {
403        for(kk=0; kk<ZDIM; kk++)
404          supertile[ii+2*IN_TILE_PTS_X][jj][kk] = newtile[ii][jj][kk];
405      }
406    }
407  }
408
409  /* Right */
410  ii = i + IN_TILE_PTS_X;
411  if (ii >= (int)(360./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X) ii -= (int)(360./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X;
412  jj = j;
413  newtile = get_tile_from_cache(ii, jj); 
414  for(ii=0; ii<IN_TILE_PTS_X; ii++)
415  {
416    for(jj=0; jj<IN_TILE_PTS_Y; jj++)
417    {
418      for(kk=0; kk<ZDIM; kk++)
419        supertile[ii+2*IN_TILE_PTS_X][jj+IN_TILE_PTS_Y][kk] = newtile[ii][jj][kk];
420    }
421  }
422
423  /* Upper-right */
424  ii = i + IN_TILE_PTS_X;
425  if (ii >= (int)(360./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X) ii -= (int)(360./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X;
426  doflip = 0;
427  jj = j + IN_TILE_PTS_Y;
428  if (jj >= (int)(180./IN_TILE_DEGREES_LAT)*IN_TILE_PTS_Y)
429  {
430    doflip = 1;
431    jj = j;
432    ii -= (int)(180./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X;
433    if (ii < 0) ii += (int)(360./IN_TILE_DEGREES_LON)*IN_TILE_PTS_X;
434  }
435  newtile = get_tile_from_cache(ii, jj); 
436  if (doflip)
437  {
438    for(ii=0; ii<IN_TILE_PTS_X; ii++)
439    {
440      for(jj=0; jj<IN_TILE_PTS_Y; jj++)
441      {
442        for(kk=0; kk<ZDIM; kk++)
443          supertile[ii+2*IN_TILE_PTS_X][jj+2*IN_TILE_PTS_Y][kk] = newtile[ii][IN_TILE_PTS_Y-1-jj][kk];
444      }
445    }
446  }
447  else
448  {
449    for(ii=0; ii<IN_TILE_PTS_X; ii++)
450    {
451      for(jj=0; jj<IN_TILE_PTS_Y; jj++)
452      {
453        for(kk=0; kk<ZDIM; kk++)
454          supertile[ii+2*IN_TILE_PTS_X][jj+2*IN_TILE_PTS_Y][kk] = newtile[ii][jj][kk];
455      }
456    }
457  }
458}
459
460int get_value(int i, int j, int k, int irad)
461{
462  int i_src, j_src;
463  int ii, jj;
464  int n, sum;
465  float r_rel;
466
467  if (i < 0)
468    i = i + IN_TILE_PTS_X*(int)(360./IN_TILE_DEGREES_LON);
469  if (i >= IN_TILE_PTS_X*(int)(360./IN_TILE_DEGREES_LON))
470    i = i - IN_TILE_PTS_X*(int)(360./IN_TILE_DEGREES_LON);
471  if (j < 0)
472    j = j + IN_TILE_PTS_Y*(int)(180./IN_TILE_DEGREES_LAT);
473  if (j >= IN_TILE_PTS_Y*(int)(180./IN_TILE_DEGREES_LAT))
474    j = j - IN_TILE_PTS_Y*(int)(180./IN_TILE_DEGREES_LAT);
475
476  i_src = i % IN_TILE_PTS_X + IN_TILE_PTS_X;
477  j_src = j % IN_TILE_PTS_Y + IN_TILE_PTS_Y;
478
479  /* Interpolate values from supertile */
480  sum = 0;
481  n = 0;
482  for(ii=i_src; ii<=i_src+irad; ii++)
483  {
484    for(jj=j_src; jj<=j_src+irad; jj++)
485    {
486      sum += supertile[ii][jj][k]; 
487      n++;
488    }
489  }
490
491  if (n > 0) return sum/n;
492  else return 0;
493}
494
495int is_in_supertile(int i, int j)
496{
497  if (i < 0)
498    i = i + IN_TILE_PTS_X*(int)(360./IN_TILE_DEGREES_LON);
499  if (i >= IN_TILE_PTS_X*(int)(360./IN_TILE_DEGREES_LON))
500    i = i - IN_TILE_PTS_X*(int)(360./IN_TILE_DEGREES_LON);
501  if (j < 0)
502    j = j + IN_TILE_PTS_Y*(int)(180./IN_TILE_DEGREES_LAT);
503  if (j >= IN_TILE_PTS_Y*(int)(180./IN_TILE_DEGREES_LAT))
504    j = j - IN_TILE_PTS_Y*(int)(180./IN_TILE_DEGREES_LAT);
505
506  /* Check whether (i,j) is in the interior of supertile */
507  if ((i >= supertile_min_x) && (i < supertile_min_x+IN_TILE_PTS_X) &&
508      (j >= supertile_min_y) && (j < supertile_min_y+IN_TILE_PTS_Y))
509    return 1;
510  else
511    return 0;
512}
513
514int main(int argc, char ** argv)
515{
516  int tile_x, tile_y, input_x, input_y, temp_x, temp_y;
517  int i, j, k, z, ii, jj;
518  int i_src, j_src;
519  int *** intdata;
520  int out_fd;
521  int ir_rel;
522  float r_rel;
523  unsigned char * outdata;
524  char out_filename[256];
525
526  r_rel = (float)IN_TILE_PTS_X/(float)OUT_TILE_PTS_X*OUT_TILE_DEGREES_LON/IN_TILE_DEGREES_LON;
527  ir_rel = (int)rint(r_rel);
528
529  /* Allocate memory to hold a single output tile */
530  intdata = (int ***)malloc(sizeof(int **)*(OUT_TILE_PTS_X+2*HALO_WIDTH)); 
531  for(i=0; i<(OUT_TILE_PTS_X+2*HALO_WIDTH); i++)
532  {
533    intdata[i] = (int **)malloc(sizeof(int *)*(OUT_TILE_PTS_Y+2*HALO_WIDTH)); 
534    for(j=0; j<(OUT_TILE_PTS_Y+2*HALO_WIDTH); j++)
535    {
536      intdata[i][j] = (int *)malloc(sizeof(int)*ZDIM); 
537    }
538  }
539
540  /* Allocate output buffer */
541  outdata = (unsigned char *)malloc((OUT_TILE_PTS_X+2*HALO_WIDTH)*(OUT_TILE_PTS_Y+2*HALO_WIDTH)*ZDIM*N_BYTES);
542
543  for(tile_x=0; tile_x<OUT_TILE_PTS_X*(int)(360.0/OUT_TILE_DEGREES_LON); tile_x+=OUT_TILE_PTS_X)
544  {
545    for(tile_y=0; tile_y<OUT_TILE_PTS_Y*(int)(180.0/OUT_TILE_DEGREES_LAT); tile_y+=OUT_TILE_PTS_Y)
546    {
547      /* Build name of output file for current tile_x and tile_y */
548      sprintf(out_filename,"retiled/%5.5i-%5.5i.%5.5i-%5.5i",tile_x+1,tile_x+OUT_TILE_PTS_X,tile_y+1,tile_y+OUT_TILE_PTS_Y);
549
550      /* Initialize the output data for current tile */
551      for(i=0; i<(OUT_TILE_PTS_X+2*HALO_WIDTH); i++) 
552      {
553        for(j=0; j<(OUT_TILE_PTS_Y+2*HALO_WIDTH); j++)
554        {
555            intdata[i][j][0] = MSG_FLAG;
556        }
557      }
558
559      /* Attempt to open output file */
560      if ((out_fd = open(out_filename, O_WRONLY|O_CREAT|O_TRUNC, S_IRUSR|S_IWUSR|S_IRGRP|S_IROTH)) == -1)
561      {
562        fprintf(stderr, "Error: Could not create or open file %s\n", out_filename);
563        return 1;
564      }
565
566      /* Fill tile with data */
567      for(i=-1*HALO_WIDTH; i<OUT_TILE_PTS_X+HALO_WIDTH; i++) 
568      {
569        for(j=-1*HALO_WIDTH; j<OUT_TILE_PTS_Y+HALO_WIDTH; j++)
570        {
571          if (intdata[i+HALO_WIDTH][j+HALO_WIDTH][0] == MSG_FLAG)
572          {
573            i_src = ir_rel*(tile_x+i); 
574            j_src = ir_rel*(tile_y+j); 
575
576            build_supertile(i_src,j_src); 
577
578            for(ii=-1*HALO_WIDTH; ii<OUT_TILE_PTS_X+HALO_WIDTH; ii++) 
579            {
580              for(jj=-1*HALO_WIDTH; jj<OUT_TILE_PTS_Y+HALO_WIDTH; jj++)
581              {
582                i_src = ir_rel*(tile_x+ii); 
583                j_src = ir_rel*(tile_y+jj); 
584                if (is_in_supertile(i_src,j_src))
585                {
586                  for(k=0; k<ZDIM; k++)
587                    intdata[ii+HALO_WIDTH][jj+HALO_WIDTH][k] = get_value(i_src,j_src,k,ir_rel-1); 
588                } 
589              }
590            }
591
592          }
593        }
594      }
595
596      /* Write out the data */
597      for(i=0; i<(OUT_TILE_PTS_X+2*HALO_WIDTH); i++) 
598      {
599        for(j=0; j<(OUT_TILE_PTS_Y+2*HALO_WIDTH); j++)
600        {
601          for(z=0; z<ZDIM; z++)
602          {
603            for(k=0; k<N_BYTES; k++)
604            {
605              outdata[z*N_BYTES*(OUT_TILE_PTS_Y+2*HALO_WIDTH)*(OUT_TILE_PTS_X+2*HALO_WIDTH)+j*N_BYTES*(OUT_TILE_PTS_X+2*HALO_WIDTH)+i*N_BYTES+k] = 
606                     (intdata[i][j][z] & (0xff << 8*(N_BYTES-k-1))) >> (8*(N_BYTES-k-1));
607            }
608          }
609        }
610      }
611      write(out_fd,(void *)outdata,(OUT_TILE_PTS_X+2*HALO_WIDTH)*(OUT_TILE_PTS_Y+2*HALO_WIDTH)*ZDIM*N_BYTES);
612      close(out_fd);
613      printf("Wrote file %s\n",out_filename);
614    }
615  }
616
617  /* Deallocate memory */
618  for(i=0; i<(OUT_TILE_PTS_X+2*HALO_WIDTH); i++)
619  {
620    for(j=0; j<(OUT_TILE_PTS_Y+2*HALO_WIDTH); j++)
621    {
622      free(intdata[i][j]);
623    }
624    free(intdata[i]);
625  }
626  free(intdata);
627
628
629  return 0;
630}
Note: See TracBrowser for help on using the repository browser.