1 /* -*- Mode: C; indent-tabs-mode: nil; c-basic-offset: 4 c-style: "K&R" -*- */
3 /*---------------------------------------------------------------------------
5 aint - Calculates mean image intensity at each interrogation area
7 Copyright (C) 2006, 2007, 2008
8 Gerber van der Graaf <gerber_graaf@users.sourceforge.net
10 This program is free software; you can redistribute it and/or modify
11 it under the terms of the GNU General Public License as published by
12 the Free Software Foundation; either version 2, or (at your option)
15 This program is distributed in the hope that it will be useful,
16 but WITHOUT ANY WARRANTY; without even the implied warranty of
17 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18 GNU General Public License for more details.
20 You should have received a copy of the GNU General Public License
21 along with this program; if not, write to the Free Software Foundation,
22 Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
24 ------------------------------------------------------------------------*/
33 #define PARFILE "gpivrc" /* Parameter file name */
35 Usage: gpiv_aint [ -h | --help] [-p | --print] [-s | --second]\n\
36 [-V | --verbose] [-v | --version] filename \n\
39 -h | --help: this on-line help \n\
40 -p | --print: print parameters to stdout \n\
41 -s | --second: obtain averages from second image frame \n\
42 -V | --verbose: prints data of interrogation area's and resulting mean \n\
43 -v | --version: version number \n\
44 filename: image file name \n\
48 gpiv_aint calculates mean image intensity at each interrogation area"
51 * Extension of filename means: "Interrogation Area Intensity"
53 #define AINT_EXT ".ain"
58 gboolean second_image
= FALSE
;
59 gboolean print_par
= FALSE
, verbose
= FALSE
;
63 * Variables for development version
71 assign_img2intarr (gint ipoint_x
,
87 assign_img2intarr_central (gint ipoint_x
,
103 int_mean (gfloat
**int_area
,
108 command_args (int argc
,
110 char fname
[GPIV_MAX_CHARS
]
114 make_fname (char *fname
,
120 static GpivScalarData
*
121 post_ia_intensity (GpivImage
*image
,
122 GpivPivData
*piv_data
,
126 static GpivScalarData
*
127 post_ia_intensity (GpivImage
*image
,
128 GpivPivData
*piv_data
,
131 /* ----------------------------------------------------------------------------
132 * Calculates mean image intensity for each I.A.
135 GpivScalarData
*sc_data
= NULL
;
140 GpivImage
*loc_image
= NULL
;
142 gboolean flag_accept
= FALSE
;
143 gint pre_shift_col_act
= 0, pre_shift_row_act
= 0;
146 assert (image
->frame1
[0] != NULL
);
147 assert (image
->frame2
[0] != NULL
);
149 assert (piv_data
->point_x
!= NULL
);
150 assert (piv_data
->point_y
!= NULL
);
151 assert (piv_data
->dx
!= NULL
);
152 assert (piv_data
->dy
!= NULL
);
153 assert (piv_data
->snr
!= NULL
);
154 assert (piv_data
->peak_no
!= NULL
);
156 int_area_1
= gpiv_matrix (piv_par
->int_size_f
, piv_par
->int_size_f
);
157 int_area_2
= gpiv_matrix (piv_par
->int_size_f
, piv_par
->int_size_f
);
160 sc_data
= gpiv_alloc_scdata (piv_data
->nx
, piv_data
->ny
);
163 * Deform entire image
165 if (piv_par
->int_scheme
== GPIV_IMG_DEFORM
) {
166 loc_image
= gpiv_cp_img (image
);
167 gpiv_imgproc_deform (loc_image
, piv_data
);
169 loc_image
->frame1
= image
->frame1
;
170 loc_image
->frame2
= image
->frame2
;
174 * For each I.A. obtain image data
175 * Source: gpiv_piv_interr_reg
177 if (piv_par
->int_scheme
== GPIV_ZERO_OFF_CENTRAL
178 || piv_par
->int_scheme
== GPIV_ZERO_OFF_FORWARD
182 g_message("post_ia_intensity:: assign_img2intarr_central");
186 pre_shift_col_act
= piv_data
->dx
[i
][j
] + piv_par
->pre_shift_col
;
187 pre_shift_row_act
= piv_data
->dy
[i
][j
] + piv_par
->pre_shift_row
;
188 for (i
= 0; i
< piv_data
->ny
; i
++) {
189 for (j
= 0; j
< piv_data
->nx
; j
++) {
191 assign_img2intarr_central ((int) piv_data
->point_x
[i
][j
],
192 (int) piv_data
->point_y
[i
][j
],
201 image
->header
->nrows
,
202 image
->header
->ncolumns
,
203 piv_par
->int_size_f
);
205 * Calculate average value and apply to out_data
207 sc_data
->point_x
[i
][j
] = piv_data
->point_x
[i
][j
];
208 sc_data
->point_y
[i
][j
] = piv_data
->point_y
[i
][j
];
210 sc_data
->flag
[i
][j
] = 1;
212 sc_data
->scalar
[i
][j
] = int_mean(int_area_2
,
213 piv_par
->int_size_f
);
215 sc_data
->scalar
[i
][j
] = int_mean(int_area_1
,
216 piv_par
->int_size_f
);
217 if (verbose
) printf("x=%d y=%d mean = %f \n\n",
218 (int)sc_data
->point_x
[i
][j
],
219 (int)sc_data
->point_x
[i
][j
],
220 sc_data
->scalar
[i
][j
]);
223 sc_data
->flag
[i
][j
] = 0;
224 sc_data
->scalar
[i
][j
] = 0.0;
231 g_message("post_ia_intensity:: assign_img2intarr");
234 pre_shift_col_act
= piv_par
->pre_shift_col
;
235 pre_shift_row_act
= piv_par
->pre_shift_row
;
236 for (i
= 0; i
< piv_data
->ny
; i
++) {
237 for (j
= 0; j
< piv_data
->nx
; j
++) {
239 assign_img2intarr((int) piv_data
->point_x
[i
][j
],
240 (int) piv_data
->point_y
[i
][j
],
249 image
->header
->nrows
,
250 image
->header
->ncolumns
,
251 piv_par
->int_size_f
);
253 * Calculate average value and apply to out_data
255 sc_data
->point_x
[i
][j
] = piv_data
->point_x
[i
][j
];
256 sc_data
->point_y
[i
][j
] = piv_data
->point_y
[i
][j
];
258 sc_data
->flag
[i
][j
] = 1;
260 sc_data
->scalar
[i
][j
] = int_mean(int_area_2
,
261 piv_par
->int_size_f
);
263 sc_data
->scalar
[i
][j
] = int_mean(int_area_1
,
264 piv_par
->int_size_f
);
265 if (verbose
) printf("x=%d y=%d mean = %f \n\n",
266 (int)sc_data
->point_x
[i
][j
],
267 (int)sc_data
->point_x
[i
][j
],
268 sc_data
->scalar
[i
][j
]);
271 sc_data
->flag
[i
][j
] = 0;
272 sc_data
->scalar
[i
][j
] = 0.0;
279 * Free allocated memory
281 gpiv_free_matrix(int_area_1
);
282 gpiv_free_matrix(int_area_2
);
283 if (piv_par
->int_scheme
== GPIV_IMG_DEFORM
) {
284 gpiv_free_img (loc_image
);
292 assign_img2intarr(gint ipoint_x
,
306 /*-----------------------------------------------------------------------------
307 * Assigns image data to the interrogation area arrays straightforward
311 gint arg_int1_r
= ipoint_y
- int_size_f
/ 2 + 1;
312 gint arg_int1_c
= ipoint_x
- int_size_f
/ 2 + 1;
313 gint arg_int2_r
= ipoint_y
- int_size_i
/ 2 + 1;
314 gint arg_int2_c
= ipoint_x
- int_size_i
/ 2 + 1;
319 g_message("assign_img2intarr:: 0 arg_int1_r = %d arg_int1_c = %d",
320 arg_int1_r
, arg_int1_c
);
321 g_message("assign_img2intarr:: 0 arg_int2_r = %d arg_int2_c = %d",
322 arg_int2_r
, arg_int2_c
);
325 * Check if Interrogation Areas are within the image boundaries.
326 * Principally arg_int1_r,c don't have to be tested as
327 * int_size_i >= int_size_f, but has been kept to maintain generallity with the
328 * other assign_img2intarr* functions
330 if ((arg_int1_r
) >= 0
331 && (arg_int1_r
+ int_size_f
- 1) < nrows
333 && (arg_int1_c
+ int_size_f
- 1) < ncolumns
336 && (arg_int2_r
+ int_size_i
- 1) < nrows
338 && (arg_int2_c
+ int_size_i
- 1) < ncolumns
) {
346 g_message("assign_img2intarr:: 1 flag_valid=%d", flag_valid
);
348 if (flag_valid
== TRUE
) {
350 * reset int_area_1, int_area_2 values
352 memset(int_area_1
[0], 0.0, (sizeof(gfloat
)) * int_size_0
* int_size_0
);
353 memset(int_area_2
[0], 0.0, (sizeof(gfloat
)) * int_size_0
* int_size_0
);
355 if (verbose
) printf("assign_img2intarr: x = %d y = %d",
357 for (m
= 0; m
< int_size_f
; m
++) {
358 if (verbose
) printf("\n");
359 for (n
= 0; n
< int_size_f
; n
++) {
361 (float) img_1
[m
+ arg_int1_r
][n
+ arg_int1_c
];
362 if (verbose
) printf("%d ", (int) int_area_1
[m
][n
]);
365 if (verbose
) printf("\n\n");
367 for (m
= 0; m
< int_size_i
; m
++) {
368 for (n
= 0; n
< int_size_i
; n
++) {
370 (float) img_2
[m
+ arg_int2_r
][n
+ arg_int2_c
];
376 g_message("assign_img2intarr:: END");
384 assign_img2intarr_central(gint ipoint_x
,
398 /*-----------------------------------------------------------------------------
399 * Assigns image data to the interrogation area arrays for central differential
404 gint idum
= gpiv_max((int_size_i
- int_size_f
) / 2, 0);
405 gint arg_int1_r
= ipoint_y
- int_size_f
/ 2 + 1 - pre_shift_row
/ 2 -
407 gint arg_int1_c
= ipoint_x
- int_size_f
/ 2 + 1 - pre_shift_col
/ 2 -
409 gint arg_int2_r
= ipoint_y
- int_size_i
/ 2 + 1 + pre_shift_row
/ 2;
410 gint arg_int2_c
= ipoint_x
- int_size_i
/ 2 + 1 + pre_shift_col
/ 2;
414 g_message("assign_img2intarr_central:: 0 arg_int1_r = %d arg_int1_c = %d",
415 arg_int1_r
, arg_int1_c
);
416 g_message("assign_img2intarr_central:: 0 arg_int2_r = %d arg_int2_c = %d",
417 arg_int2_r
, arg_int2_c
);
420 * Check if Interrogation Areas are within the image boundaries.
422 if ((arg_int1_r
) >= 0
423 && (arg_int1_r
+ int_size_f
- 1) < nrows
425 && (arg_int1_c
+ int_size_f
- 1) < ncolumns
428 && (arg_int2_r
+ int_size_i
- 1) < nrows
430 && (arg_int2_c
+ int_size_i
- 1) < ncolumns
) {
438 g_message("assign_img2intarr_central:: 1 flag_valid=%d", flag_valid
);
440 if (flag_valid
== TRUE
) {
442 * reset int_area_1, int_area_2 values
444 memset(int_area_1
[0], 0.0, (sizeof(gfloat
)) * int_size_0
* int_size_0
);
445 memset(int_area_2
[0], 0.0, (sizeof(gfloat
)) * int_size_0
* int_size_0
);
447 if (verbose
) printf("assign_img2intarr_central: x = %d y = %d",
449 for (m
= 0; m
< int_size_f
; m
++) {
450 if (verbose
) printf("\n");
451 for (n
= 0; n
< int_size_f
; n
++) {
452 int_area_1
[m
+ idum
][n
+ idum
] =
453 (float) img_1
[m
+ arg_int1_r
][n
+ arg_int1_c
];
454 if (verbose
) printf(" %d ", (int) int_area_1
[m
][n
]);
457 if (verbose
) printf("\n\n");
459 for (m
= 0; m
< int_size_i
; m
++) {
460 for (n
= 0; n
< int_size_i
; n
++) {
462 (float) img_2
[m
+ arg_int2_r
][n
+ arg_int2_c
];
469 g_message("assign_img2intarr_central:: END");
477 int_mean(gfloat
**int_area
,
480 /* ----------------------------------------------------------------------------
481 * calculates mean value from interrogation area
485 gfloat int_area_sum
= 0;
488 for (m
= 0; m
< int_size
; m
++) {
489 for (n
= 0; n
< int_size
; n
++) {
490 int_area_sum
+= int_area
[m
][n
];
493 mean
= int_area_sum
/ (int_size
* int_size
);
501 command_args (int argc
,
503 char fname
[GPIV_MAX_CHARS
]
505 /* ----------------------------------------------------------------------------
506 * Command line argument handling
513 while (--argc
> 0 && (*++argv
)[0] == '-') {
516 * argc_next is set to 1 if the next cmd line argument has to be searched for;
517 * in case that the command line argument concerns more than one char or cmd
518 * line argument needs a parameter
520 while (argc_next
== 0 && (c
= *++argv
[0])) {
525 * Use Git revision control system
528 printf ("git hash: %s\n", GIT_REV
);
530 printf ("version: %s\n", GPIVTOOLS_VERSION
);
536 printf ("%s\n", argv
[0]);
537 printf ("%s\n", HELP
);
538 printf ("%s\n", USAGE
);
558 if (strcmp ("-help", *argv
) == 0) {
559 printf ("\n%s", argv
[0]);
560 printf ("\n%s", HELP
);
561 printf ("\n%s", USAGE
);
563 } else if (strcmp ("-print", *argv
) == 0) {
565 } else if (strcmp ("-second", *argv
) == 0) {
567 } else if (strcmp ("-verbose", *argv
) == 0) {
569 } else if (strcmp ("-version", *argv
) == 0) {
571 printf ("git hash: %s\n", GIT_REV
);
573 printf ("version: %s\n", GPIVTOOLS_VERSION
);
577 gpiv_error ("%s: unknown option: %s", argv
[0], *argv
);
583 fprintf (stderr
, USAGE
);
592 strcpy (fname
, *argv
);
594 gpiv_error ("%s: %s", argv
[0], USAGE
);
602 make_fname (char *fname
,
607 /* ---------------------------------------------------------------------------
608 * function to generate filenames
611 gchar
*err_msg
= NULL
;
612 gchar
*fname_base
= NULL
;
614 if (fname
== NULL
) {
615 err_msg
= "make_fname: \"fname == NULL\"";
618 fname_base
= g_strdup (fname
);
619 strtok (fname_base
, ".");
621 gpiv_io_make_fname (fname_base
, GPIV_EXT_PAR
, fname_par
);
622 if (print_par
) printf("# Parameter file: %s\n",fname_par
);
624 gpiv_io_make_fname (fname_base
, GPIV_EXT_PIV
, fname_piv
);
625 if (print_par
) printf("# PIV data file: %s\n",fname_piv
);
627 gpiv_io_make_fname (fname_base
, AINT_EXT
, fname_out
);
628 if (print_par
) printf("# Output file: %s\n",fname_out
);
640 /* ----------------------------------------------------------------------------
641 * Start of the main program
644 gchar
* err_msg
= NULL
;
646 gchar fname
[GPIV_MAX_CHARS
], fname_out
[GPIV_MAX_CHARS
],
647 fname_par
[GPIV_MAX_CHARS
], fname_piv
[GPIV_MAX_CHARS
];
649 GpivPivPar
*piv_par
= g_new0 (GpivPivPar
, 1);
650 GpivPostPar
*post_par
= g_new0 (GpivPostPar
, 1);
652 GpivImage
*image
= NULL
;
653 GpivPivData
*piv_in_data
= NULL
;
654 GpivScalarData
*sc_out_data
= NULL
;
658 gpiv_piv_parameters_set (piv_par
, FALSE
);
659 gpiv_post_parameters_set (post_par
, FALSE
);
662 command_args (argc
, argv
, fname
);
666 * Generating proper filenames
669 make_fname (fname
, fname_par
, fname_piv
, fname_out
))
671 gpiv_error ("%s: Failure calling make_fname", argv
[0]);
676 * Prints command line parameters to par-file
678 if ((fp_par
= fopen (fname_par
, "a")) == NULL
) {
679 gpiv_error ("\n%s: failure opening %s for input",
683 fprintf (fp_par
, "\n\n# Software: %s\n# git hash: %s\n# Command line options:\n",
686 fprintf (fp_par
, "\n\n# Software: %s\n# version: %s\n# Command line options:\n",
687 argv
[0], GPIVTOOLS_VERSION
);
691 * Reading parametes from PARFILE and (system) resources (and writing to data
695 gpiv_scan_parameter (GPIV_PIVPAR_KEY
, PARFILE
, piv_par
, print_par
);
697 gpiv_scan_resourcefiles (GPIV_PIVPAR_KEY
, piv_par
, print_par
))
698 != NULL
) gpiv_error ("%s: %s", argv
[0], err_msg
);
699 gpiv_piv_print_parameters (fp_par
, piv_par
);
701 gpiv_scan_parameter (GPIV_POSTPAR_KEY
, PARFILE
, post_par
, print_par
);
703 gpiv_scan_resourcefiles (GPIV_POSTPAR_KEY
, post_par
, print_par
))
704 != NULL
) gpiv_error ("%s: %s", argv
[0], err_msg
);
705 gpiv_post_print_parameters (fp_par
, post_par
);
709 gpiv_post_check_parameters_read (post_par
, NULL
))
710 != NULL
) gpiv_error ("%s: %s", argv
[0], err_msg
);
715 if ((image
= gpiv_fread_image (fname
)) == NULL
) {
716 gpiv_error ("%s: gpiv_fread_image_ni", argv
[0]);
721 * Read PIV data for input
723 if ((fp
= fopen (fname_piv
, "rb")) == NULL
) {
726 fprintf (fp_par
, "# WARNING: no PIV data found; particle displacements set to 0");
727 if (print_par
) gpiv_warning ("Particle displacements set to 0");
730 gpiv_piv_count_pivdata_fromimage (image
->header
, piv_par
, &nx
, &ny
))
731 != NULL
) gpiv_error ("%s: %s", argv
[0], err_msg
);
734 gpiv_piv_gridgen (nx
, ny
, image
->header
, piv_par
))
735 == NULL
) gpiv_error ("%: %s", argv
[0], err_msg
);
738 gpiv_0_pivdata (piv_in_data
))
739 != NULL
) gpiv_error (": %s", argv
[0], err_msg
);
742 if ((piv_in_data
= gpiv_read_pivdata (fp
))
743 == NULL
) gpiv_error ("%s: %s", argv
[0], err_msg
);
751 * Here the function calls of the post-processing; intensity of
752 * Interrogation Area's
754 if (print_par
== TRUE
) printf ("\n");
757 post_ia_intensity (image
, piv_in_data
, piv_par
))
758 == NULL
) gpiv_error ("%s: Failure %s", argv
[0], err_msg
);
761 * And writing to output
763 if ((fp
= fopen (fname_out
, "wb")) == NULL
) {
764 gpiv_error ("%s: failing opening %s", argv
[0], fname_out
);
767 gpiv_write_scdata (fp
, sc_out_data
, TRUE
))
768 != NULL
) gpiv_error ("%s: %s", argv
[0], err_msg
);
772 * Freeing allocated memory
774 gpiv_free_pivdata (piv_in_data
);
775 gpiv_free_img (image
);
776 if (print_par
== TRUE
) printf ("\n");