For most of the naturally occurring signals, especially biomedical signals, the underlying physical process generating the signal is often not fully known, making it difficult to obtain a parametric model. Therefore, signal processing techniques are used to analyze the signal for non-parametrically characterizing the underlying system from which the signals are produced. Most of the real life systems are nonlinear and time varying, which poses a challenge while characterizing them. Additionally, multiple sensors are used to extract signals from such systems, resulting in multichannel signals which are inherently coupled. In this thesis, we counter this challenge by using Recurrence Plot based techniques for characterizing biomedical systems such as heart or brain, using signals such as heart rate variability (HRV), electroencephalogram(EEG) or functional magnetic resonance imaging (fMRI), respectively, extracted from them. In time series analysis, it is well known that a system can be represented by a trajectory in an N-dimensional state space, which completely represents an instance of the system behavior. Such a system characterization has been done using dynamical invariants such as correlation dimension, Lyapunov exponent etc. Takens has shown that when the state variables of the underlying system are not known, one can obtain a trajectory in ‘phase space’ using only the signals obtained from such a system. The phase space trajectory is topologically equivalent to the state space trajectory. This enables us to characterize the system behavior from only the signals sensed from them. However, estimation of correlation dimension, Lyapunov exponent, etc, are vulnerable to non-stationarities in the signal and require large number of sample points for accurate computation, both of which are important in the case of biomedical signals. Alternatively, a technique called Recurrence Plots (RP) has been proposed, which addresses these concerns, apart from providing additional insights. Measures to characterize RPs of single and two channel data are called Recurrence Quantification Analysis (RQA) and cross RQA (CRQA), respectively. These methods have been applied with a good measure of success in diverse areas. However, they have not been studied extensively in the context of experimental biomedical signals, especially multichannel data. In this thesis, the RP technique and its associated measures are briefly reviewed. Using the computational tools developed for this thesis, RP technique has been applied on select single channel, multichannel and multimodal (i.e. multiple channels derived from different modalities) biomedical signals. Connectivity analysis is demonstrated as post-processing of RP analysis on multichannel signals such as EEG and fMRI. Finally, a novel metric, based on the modification of a CRQA measure is proposed, which shows improved results. For the case of single channel signal, we have considered a large database of HRV signals of 112 subjects recorded for both normal and abnormal (anxiety disorder and depression disorder) subjects, in both supine and standing positions. Existing RQA measures, Recurrence Rate and Determinism, were used to distinguish between normal and abnormal subjects with an accuracy of 58.93%. A new measure, MLV has been introduced, using which a classification accuracy of 98.2% is obtained. Correlation between probabilities of recurrence (CPR) is a CRQA measure used to characterize phase synchronization between two signals. In this work, we demonstrate its utility with application to multimodal and multichannel biomedical signals. First, for the multimodal case, we have computed running CPR (rCPR), a modification proposed by us, which allows dynamic estimation of CPR as a function of time, on multimodal cardiac signals (electrocardiogram and arterial blood pressure) and demonstrated that the method can clearly detect abnormalities (premature ventricular contractions); this has potential applications in cardiac care such as assisted automated diagnosis. Second, for the multichannel case, we have used 16 channel EEG signals recorded under various physiological states such as (i) global epileptic seizure and pre-seizure and (ii) focal epilepsy. CPR was computed pair-wise between the channels and a CPR matrix of all pairs was formed. Contour plot of the CPR matrix was obtained to illustrate synchronization. Statistical analysis of CPR matrix for 16 subjects of global epilepsy showed clear differences between pre-seizure and seizure conditions, and a linear discriminant classifier was used in distinguishing between the two conditions with 100% accuracy. Connectivity analysis of multichannel EEG signals was performed by post-processing of the CPR matrix to understand global network-level characterization of the brain. Brain connectivity using thresholded CPR matrix of multichannel EEG signals showed clear differences in the number and pattern of connections in brain connectivity graph between epileptic seizure and pre-seizure. Corresponding brain headmaps provide meaningful insights about synchronization in the brain in those states. K-means clustering of connectivity parameters of CPR and linear correlation obtained from global epileptic seizure and pre-seizure showed significantly larger cluster centroid distances for CPR as opposed to linear correlation, thereby demonstrating the efficacy of CPR. The headmap in the case of focal epilepsy clearly enables us to identify the focus of the epilepsy which provides certain diagnostic value. Connectivity analysis on multichannel fMRI signals was performed using CPR matrix and graph theoretic analysis. Adjacency matrix was obtained from CPR matrices after thresholding it using statistical significance tests. Graph theoretic analysis based on communicability was performed to obtain community structures for awake resting and anesthetic sedation states. Concurrent behavioral data showed memory impairment due to anesthesia. Given the fact that previous studies have implicated the hippocampus in memory function, the CPR results showing the hippocampus within the community in awake state and out of it in anesthesia state, demonstrated the biological plausibility of the CPR results. On the other hand, results from linear correlation were less biologically plausible. In biological systems, highly synchronized and desynchronized systems are of interest rather than moderately synchronized ones. However, CPR is approximately a monotonic function of synchronization and hence can assume values which indicate moderate synchronization. In order to emphasize high synchronization/ desynchronization and de-emphasize moderate synchronization, a new method of Correlation Synchronization Convergence Time (CSCT) is proposed. It is obtained using an iterative procedure involving the evaluation of CPR for successive autocorrelations until CPR converges to a chosen threshold. CSCT was evaluated for 16 channel EEG data and corresponding contour plots and histograms were obtained, which shows better discrimination between synchronized and asynchronized states compared to the conventional CPR. This thesis has demonstrated the efficacy of RP technique and associated measures in characterizing various classes of biomedical signals. The results obtained are corroborated by well known physiological facts, and they provide physiologically meaningful insights into the functioning of the underlying biological systems, with potential diagnostic value in healthcare.