\section{Introduction}
\label{sec:intro}
Millisecond pulsars (MSPs) are fast rotating neutron stars with exceptional rotational stability enabling the precise determination of their rotational and orbital (for systems in binary) properties as well as using them to probe the interstellar medium (ISM). Pulsar timing studies enable modeling of the deviations of observed pulse arrival times from the values predicted by timing models while pulsar signal propagates through dynamic environment of the ISM. The stochastic gravitational wave (GW) background manifests as an un-modeled effect in the timing residuals whose detectability depends on the timing span and precision of the measurements \citep{2013CQGra..30v4015S}.  The Pulsar Timing Array (PTA) experiment \citep{2021arXiv210110081V} uses a set of MSPs with different angular separations in the sky to search for the angular correlation between the residuals of the arrival times of pairs of pulsars (i.e., Hellings and Downs curve, \cite{2015AmJPh..83..635J}). Such correlation reveals the signature of low-frequency stochastic isotropic GW background in the timing data, where the largest contribution is thought to be coming from an ensemble of merging super-massive black hole binaries \citep{2019A&ARv..27....5B}. 

One of the crucial challenges for PTAs is to disentangle and mitigate the timing noise contributed by variations in the free electron density of the interstellar plasma. Time-varying ISM effects (i.e., changes in the dispersion measure, the influence of scattering) on pulse arrival time need to be precisely determined to improve the timing precision. 
The emission from the pulsar undergoes frequency dependent effects as it propagates through the ISM. A signal of frequency $\nu$ arrives at Earth at a delayed time $\Delta t_{\nu}$, with respect to infinite frequency, which is given by
\begin{equation}
    \Delta t_{\nu}=K\times DM\,\nu^{-2}
\end{equation}
where K is a constant with value of 4.148808(3) GHz$^2\,$cm$^3\,$pc$^{-1}\,$ms and DM (dispersion measure) is the free electron column density along the line-of-sight (LOS)
\begin{equation}
    DM\equiv \int_{LOS} n_e dl  
\end{equation}
 where the integration is along the LOS from observatory to pulsar. Equation (1) shows that a typical DM variation of $10^{-3}-10^{-4}$\,pc  cm$^{-3}$, seen in pulsar observations \citep{2020A&A...644A.153D},
 introduces a change in pulse arrival time of more than $1 \mu s$ at $\nu\sim$ 1 GHz (with respect to the infinite frequency). Whereas, to achieve timing accuracy better than 100 ns at an observing frequency of 1400 MHz, the DM variation needs to be modeled at a precision of $\sim$ $10^{-5}$\,pc cm$^{-3}\,$ \citep{2007MNRAS.378..493Y}. 
 
Data collected by the International Pulsar Timing Array (IPTA) consists of several MSPs ($\sim65$, \cite{2019MNRAS.490.4666P}) observed over a wide range of frequencies ($0.3-3.1$ GHz) from various telescopes with the aim to improve the
limits on the isotropic stochastic GW background. Due to the weakening of ISM effects, frequencies greater than 1 GHz are preferred for high-precision timing analysis. However, low-frequencies (i.e., $<$ 1 GHz) can provide a sensitive probe for measuring DM and its temporal evolution to mitigate adverse effects in the arrival times (\cite{2012A&A...543A..66H}) which are embedded in the high frequency measurements.

The upgraded GMRT (uGMRT, \cite{2017CSci..113..707G}) is one of the most sensitive radio telescopes at low radio frequencies and it can cover a frequency range from 120 to 1460 MHz. The GMRT being a IPTA telescope can provide sensitive low-frequency timing measurements, which are already demonstrated by \cite{2020arXiv200908409J} with the legacy GMRT and \cite{2021arXiv210105334K} with the uGMRT. Band-3 of the uGMRT (i.e. 300 to 500 MHz) with its large fractional bandwidth provides a facility for very accurate in-band DM estimates. Assuming frequency independent DM \citep{2016ApJ...817...16C}, precise DM measurements obtained from this band can be used to correct for dispersive delays in simultaneous high-frequency timing data. According to the radiometer equation \citep{2012hpa..book.....L}, a larger observing bandwidth results in a higher signal-to-noise ratio (S/N) pulse profile promising better DM precision. However, at low frequencies, radio frequency interference (RFI), scintillation, and scattering affect the pulse detection significance. Also, the intrinsic pulse profile can evolve significantly with frequency. The standard narrow-band timing technique \citep{2012hpa..book.....L} uses a single frequency averaged template to generate times-of-arrival (ToAs) for different sub-bands within the observing bandwidth. The narrow-band technique doesn't account for any frequency-dependent effects. It works adequately at high frequencies with smaller fractional bandwidth where the frequency-dependent effects are negligible. \cite{2014ApJ...790...93P} (along with \cite{2014MNRAS.443.3752L}) describes the simultaneous wide-band ToA and DM measurement technique using a frequency-dependent template.  \cite{2019ApJ...871...34P} developed a principal-component-decomposition based modeling of pulse profiles as a function of frequency, which is an input to the wide-band ToA+DM measurement technique. All of these are implemented in a package called ``$\it{PulsePortraiture}$"\footnote{\label{pulse_portraiture}\url{https://github.com/pennucci/PulsePortraiture}} \citep{2016ascl.soft06013P}. Using {\textit{PulsePortraiture}$^{\ref{pulse_portraiture}}$}, we can estimate the ToA and DM simultaneously at a high precision with frequency-dependent template.

As it's first real application, \cite{2021ApJS..252....5A} reported wide-band timing results for 47 NANOGrav MSPs with a range of frequency coverage: 1.4 GHz (with a bandwidth of $\sim$600 MHz), 800 MHz (with a bandwidth of $\sim$ 186 MHz) and 430 MHz (with a bandwidth 25-50 MHz). The detailed comparisons with narrow-band timing results for these MSPs establish the potential of wide-band timing in achieving higher timing and DM precision. 
In this work, we present the results of applying wide-band timing analysis for four GMRT discovered MSPs (having flux densities around a mJy) with the uGMRT in band-3 (300-500 MHz) and band-4 (550-750 MHz). We also highlight results for a few bright PTA MSPs observed with the uGMRT in band-3 and band-5 (1060-1460 MHz) for validation of the wide-band timing pipeline. Since the frequency-dependent effects are much more prominent at low frequencies, the observing bands of uGMRT, specially band-3 with 0.5 fractional bandwidth, demonstrate effectiveness of pulse-profile modeling with frequency. Observation and data processing details are provided in Section \ref{sec:post}. Section \ref{sec:techniques} contains the details of the narrow-band and wide-band timing techniques. Section \ref{sec:Results} contains the measurements obtained from two timing analysis and compare that with some of the existing results (\cite{2021arXiv210105334K}, \cite{2021ApJS..252....5A} and \cite{2021arXiv211206908N}). In Section \ref{sec:conclusions} we summarise the improvements seen with the wide-band timing analysis.

\section{Observations and data processing}
\label{sec:post}

We observed 4 GMRT discovered pulsars (from now on we will refer them as ``non-PTA pulsars") and 4 PTA pulsars, using an instantaneous bandwidth of 200 MHz for band-3 (300$-$500 MHz), band-4 (550$-$750 MHz), and 400 MHz bandwidth for band-5 (1060$-$1460 MHz). Table \ref{table_flux} lists the period, DM, and flux densities of these eight MSPs.

The phased array beam of the uGMRT was recorded after online coherent dedispersion (where each sub-band voltage samples are corrected for dispersive delays due to the ISM), for band-3 and band-4, with a spectral resolution of 390 kHz at 10.24 $\mu$s time-resolution. In parallel, Stokes-I filterbank data at 48 kHz spectral and 81.92 $\mu$s time-resolution was also acquired for offline incoherent dedispersion (where the intra-channel dispersion smearing is not corrected). In band-3, we mask the 360$-$380 MHz frequency band affected by the persistent MUOS emission. Table \ref{table_1} shows the S/N of the observed MSPs in different frequency bands of the uGMRT. The non-PTA pulsars have good S/Ns in band-3. Among them, PSR J1646$-$2142 and J2144$-$5237 are also bright in band-4. 
PSR J1120$-$3618 and J1828$+$0625 are having relatively lower detection significance in band-4, thus excluded from this work. All of the non-PTA pulsars show significant profile evolution with frequency (Fig. \ref{GMRT_pulsar_profiles}) and therefore provide a wide-variety sample to check the efficiency  of wide-band timing. PTA pulsars J2145$-$0750, J1640$+$2224, J1713$+$0747, and J1909$-$3744 are some of the best-timed MSPs and were chosen based on the recent high precision timing results with the legacy GMRT reported in \cite{2020arXiv200908409J}. Long-term timing studies of these pulsars are available in \cite{2021ApJS..252....5A} allowing us to evaluate wide-band timing with the uGMRT at low frequency. We have also included band-5 observations for these MSPs to compare with the corresponding timing precision achieved in \cite{2021ApJS..252....5A}.
\begin{table}[H]
\begin{center}
\begin{tabular}{|c|c|c|c c c|}
 \hline
  &&Dispersion&\multicolumn{3}{c|}{Flux density $^{\dagger}$} \\*
  PSR&Period&Measure& 400 MHz & 650 MHz & 1260 MHz \\*
  &(ms)&($cm^{-3}\,pc$)& \multicolumn{3}{c|}{(mJy)}\\*
  \hline\hline
  {J1646$-$2142} & 5.85 & 29.74 & 2.2 & 1.1 & -  \\*
  {J2144$-$5237} & 5.04 & 19.55 & 1.2 & 0.6 & - \\*
  {J1120$-$3618} & 5.56 & 45.13 & 1.1 & - & -  \\*
  {J1828$+$0625} & 3.63 & 22.42 &1.3 & - & -  \\*
  \hline
  J2145$-$0750 &16.05 &  9.00 & 46$^{*}$ &-& 12.6 \\*
  J1640$+$2224 & 3.16 & 18.43 & 28.1 &-& 2.5 \\*
  J1713$+$0747 & 4.57 & 15.98 & 6.8$^{*}$ &-& 10.3 \\*
  J1909$-$3744 & 2.95 & 10.39 & 10.1 &-& 2.9 \\*
  \hline
  \end{tabular}
  \caption{The table lists the basic parameters of the observed pulsars. The flux density values and spectral indicies for the first four pulsars (GMRT discovered) are from \cite{2019ApJ...881...59B} and Bhattacharyya et al. (in prep). For the rest four PTA pulsars, the values are available in ATNF Pulsar Catalogue\protect\footnotemark $\,$ \citep{2005AJ....129.1993M}.\\*
  $\dagger\quad$ Assuming a power-law spectrum, we estimate the flux density values at 400, 650, and 1260 MHz from the spectral index and the available flux values at nearby frequencies.\\*
  $*\quad$ Available in ATNF pulsar catalogue$^{\ref{footnotemark1}}$.}
  \label{table_flux}
  \end{center}
  \end{table}
  \vspace{-10mm}
  The four PTA pulsars were observed once a month and covered a span of $\sim$ nine months, while more than 2 years of data were available for all non-PTA pulsars. Fig. \ref{fig:1} shows the cadence for all observed pulsars in the different frequency bands. Since the aim was to observe with the maximum time-domain sensitivity of the uGMRT, all of these were single sub-array observations, where $\sim$ 70\% of the GMRT array was phased and combined to form a single dish with an equivalent gain of $\sim$ 7 K/Jy. A bandwidth of 200 MHz, from 300-500 MHz and 550-750 MHz, was used for band-3 and band-4 observations, respectively. Our observations in band-5 have a frequency range from 1060 to 1460 MHz. We performed online coherent de-dispersion at the best-known DM values in band-3 and band-4, whereas band-5 data was incoherently-dispersed (as online coherent dedispersion mode is available up to 200 MHz bandwidth only).
  \footnotetext{\label{footnotemark1}\url{https://www.atnf.csiro.au/research/pulsar/psrcat/}}\\*
  We also included incoherently de-dispersed data available for MSPs J1646$-$2142, J2144$-$5237, and J1120$-$3618 in band-3 and band-4 taken before online coherent de-dispersion mode was established.

\begin{table}[H]
 \begin{tabular}{||c|c|c|c|c||}
 \hline
 MSP & Mean Observation & Median S/N & No. of & Timing base-\\ [0.5ex]
 & time (mins) & Band 3~~~Band 4~~~Band 5 & Epochs &-line (years)\\
 \hline\hline
 J1646$-$2142 & 40 &  80~~~~~~~~40~~~~~~~~-- & 39 &4.08\\
 J2144$-$5237 & 55 &  50~~~~~~~~80~~~~~~~~-- & 28 &4.00\\
 J1120$-$3618 & 50 & 70~~~~~~~~~-~~~~~~~~~-- & 19 &3.25\\
 J1828$+$0625 & 40 & 30~~~~~~~~~-~~~~~~~~~-- & 15 &2.08\\
 \hline
 J2145$-$0750 & 30 &  2600~~~~~~~~-~~~~~~~~200 & 15 &0.75\\
 J1640$+$2224 & 25 &  400~~~~~~~~~-~~~~~~~~~60 & 12 &0.75\\
 J1713$+$0747 & 25 &  160~~~~~~~~~-~~~~~~~~250 & 12 &0.75\\
 J1909$-$3744 & 20 &  120~~~~~~~~~-~~~~~~~~~20 & 7 & 0.75\\
 \hline
\end{tabular}
\caption{The observational properties for the GMRT discovered pulsars (first four) and PTA pulsars (last four) used in this work, along with their detection significance.}
\label{table_1}
\end{table}

\begin{figure}[H]
    \centering
    \includegraphics[width=16cm, keepaspectratio]{Cadence_of_observation.png}
    \caption{The cadence and length of the observing campaign for each pulsar in different uGMRT bands (marked by separate colours) for different observing modes, where C stands for observations with online coherent dedispersion and I for offline incoherent dedispersion (following Fig. 1 of \cite{2021ApJS..252....5A}).}
    \label{fig:1}
\end{figure}


The uGMRT wide-band \citep{2017JAI.....641011R} filterbank data was de-dispersed and folded using PRESTO \citep{2011ascl.soft07017R}. The output folded-data-cubes were converted to FITS format for further analysis using PSRCHIVE \citep{2012AR&T....9..237V}. For DM estimation, we divided individual uGMRT bands into 16 sub-bands, with 128/64 bins over pulse longitudes.\\*
Fig. \ref{PTA_pulsar_profiles} shows the profile evolution with frequency and pulse profiles at 4 sub-bands for the PTA pulsars in band-3 and band-5 of the uGMRT. Fig. \ref{GMRT_pulsar_profiles} shows the profile evolution with frequency using sub-band pulse profiles of the non-PTA pulsars in band-3 and band-4. 

\begin{figure}[H]
    \subfloat[J2145-0750\label{a}]{{\includegraphics[width=0.243\linewidth,height=5cm]{J2145-0750_portrait_band3.png} }}
    \subfloat[J1640+2224\label{b}]{{\includegraphics[width=0.243\linewidth,height=5cm]{J1640+2224_portrait_band3.png} }}
    \subfloat[J1909-3744\label{c}]{{\includegraphics[width=0.243\linewidth,height=5cm]{J1909-3744_portrait_band3.png} }}
    \subfloat[J1713+0747\label{d}]{{\includegraphics[width=0.243\linewidth,height=5cm]{J1713-0747_portrait_band3.png} }}



    \subfloat[J2145-0750\label{e}]{{\includegraphics[width=0.243\linewidth,height=5cm]{J2145-0750_portrait_band5.png} }}
    \subfloat[J1640+2224\label{f}]{{\includegraphics[width=0.243\linewidth,height=5cm]{J1640+2224_portrait_band5.png} }}
    \subfloat[J1909-3744\label{g}]{{\includegraphics[width=0.243\linewidth,height=5cm]{J1909-3744_portrait_band5.png} }}
    \subfloat[J1713+0747\label{h}]{{\includegraphics[width=0.243\linewidth,height=5cm]{J1713-0747_portrait_band5.png} }}
    \caption{Pulse profiles of the PTA pulsars in band-3 and band-5 of the uGMRT showing evidences of significant profile evolution. We combined all observations in each of the uGMRT band separately to create these template profiles. The profile traces of 4 sub-bands for each of the pulsars are added to highlight the profile evolution.}
    \label{PTA_pulsar_profiles}
    \end{figure}
    \begin{figure}[H]
    \subfloat[J1646-2142\label{j}]{{\includegraphics[width=0.243\linewidth,height=5cm]{J1646-2142_portrait_band3.png} }}
    \subfloat[J2144-5237\label{k}]{{\includegraphics[width=0.243\linewidth,height=5cm]{J2144-5237_portrait_band3.png} }}
    \subfloat[J1120-3618\label{l}]{{\includegraphics[width=0.243\linewidth,height=5cm]{J1120-3618_portrait_band3.png} }}
    \subfloat[J1828+0625\label{i}]{{\includegraphics[width=0.243\linewidth,height=5cm]{J1828+0625_portrait_band3.png} }}
    
    
    \subfloat[J1646-2142\label{m}]{{\includegraphics[width=0.243\linewidth,height=5cm]{J1646-2142_portrait_band4.png} }}
    \subfloat[J2144-5237\label{n}]{{\includegraphics[width=0.243\linewidth,height=5cm]{J2144-5237_portrait_band4.jpg} }}
    \caption{Pulse profiles for  non-PTA pulsars (similar to PTA pulsars shown in Fig. 2) in band-3 and band-4 of the uGMRT. Significant profile evolution with frequency and signatures of scintillation are seen in the uGMRT band. Since PSR J1120$-$3618 and J1828$+$0625 are having lower detection significance in band-4, we haven't used that for the current work.}
    \label{GMRT_pulsar_profiles}
\end{figure}

\section{Timing techniques}
\label{sec:techniques}
\subsection{Narrow-band timing}
In narrow-band (NB) timing analysis, all the observations, with significant pulse detection, are aligned using frequency invariant phase offsets. The aligned data is then averaged in time and frequency to create a reference template profile. We use different templates for individual frequency bands of the uGMRT. Gaussians are fitted, using PSRCHIVE \citep{2012AR&T....9..237V}, to these averaged profiles to create analytic noise-free templates. We average out each epoch data in time while keeping the frequency resolution (divided into 16 frequency sub-bands for all bands) and determine pulse arrival times for different frequency sub-bands.
All narrow-band ToAs are estimated using the technique described in Appendix A of \cite{1992PTRSL.341..117T}. For a given band, we use a single profile template to determine the ToAs for each of the different frequency sub-bands.

ToA generation and model parameter fits are two separate steps in NB analysis. For measuring the temporal variation of DM, it is fitted individually for each epoch, keeping other parameters fixed in TEMPO2 \citep{2012ascl.soft10015H}. 

\subsection{Wide-band timing}

The NB technique can give accurate results when the fractional bandwidth is small or pulse shape variation with frequency is minimal. However, for the large bandwidths at low frequency (i.e., large fractional bandwidths), one needs to account for profile evolution with frequency. 
We use the required notation as described in \cite{2014ApJ...790...93P}, \cite{2019ApJ...871...34P} and include them here for completeness. 
Assuming a time-domain model has the form 
\begin{equation}
    D(\nu, \phi)=B(\nu)+a(\nu)P(\nu, \phi-\phi(\nu))+N(\nu, \phi),
\end{equation}
where $a(\nu)$ and $\phi(\nu)$ are the required scale and shift of the data $D(\nu, \phi)$ from a template $P(\nu, \phi-\phi(\nu))$. $B(\nu)$ is the effective band-pass shape of the receiver and $N(\nu, \phi)$ is additive noise assumed to be normally distributed in the absence of RFI.


The receiver band-width is divided into $n$ frequency channels with the nth channel having central frequency $\nu_n$.
A pulse profile at frequency $\nu_n$ is sampled into equally spaced intervals in pulse-phase. The one-dimensional Discrete Fourier Transform (DFT)  of equation (3), after discretising the terms, with respect to rotational phase $\phi$, and making use of the discrete Fourier shift theorem leads to
\begin{equation}
    d_{nk}=a_n p_{nk}e^{-2\pi i k \phi_n} + n_{nk},
\end{equation}
where k is the Fourier frequency, n is nth frequency channel. $d_{nk}$, $p_{nk}$ and $n_{nk}$ are the DFT of the data D, the template P and the noise N defined in equation (3). 
Minimising 
        \begin{equation}
             \chi^2(\phi_n, a_n)=\sum_{n,k}\frac{|d_{nk}-a_nt_{nk}e^{-2 i \pi k \phi_n}|^2}{\sigma^2_n}
             \label{chi_sq_eqn_WB}
        \end{equation}
        will provide estimates of the scaling factor $a_n$ and phase shift $\phi_n$ between the data $d_{nk}$ and the template $t_{nk}$. Each term in the sum is weighted by the square of the noise $\sigma_n$ estimated in the corresponding frequency channel. This approach is similar to the technique described in \cite{1992PTRSL.341..117T}, but in this case, the template can evolve with frequency. The wide-band (WB) technique is enabled to simultaneously measure the ToA ($\phi_o$ in Fourier domain) and DM by inclusion of the constraint in the equation below
        \begin{equation}
            \phi_n(\nu_n)=\phi_o+\frac{K\times DM}{P_s}(\nu_n^{-2}-\nu_{\phi_o}^{-2})
            \label{equation_fit}
        \end{equation} 
        where $\phi_o$ is the constant phase offset independent of DM, $P_s$ is the period of the pulsar, K is the dispersion constant, and $\nu_{\phi_o}$ is a choice of parameterization. The software gives freedom to select the value of $\nu_{\phi_o}$. However, we have used the default feature of the package which finds the value of $\nu_{\phi_o}$ such that there is zero covariance between $\phi_{o}$ and DM. \cite{2014ApJ...790...93P} gives a more comprehensive description of the wide-band ToA generation technique.
        
        
        Here we provide a brief description of the template creation procedure in wide-band timing, and refer to \cite{2019ApJ...871...34P} for further details. In our analysis, we select folded archives having high S/N, with same central frequency to create the profile template. Bad frequency channels, in each folded data cube, are zapped manually. We sum all the data sets in time and create 16 sub-bands over a bandwidth of 200 MHz in band-3, 4 and 400 MHz in band-5. All the data sets for a given frequency band, and pulsar, are aligned by determining a constant offset + offset proportional to $\nu^{-2}$ relative to an initial guess, where the $\nu^{-2}$ factor accounts for DM variability from one observation to another. The aligned data is averaged together. Upon iteration, it uses that result as the new reference for alignment and the process is repeated multiple times to create an ``average portrait''. Note that, here we maintain the frequency resolution during alignment. Then the average portrait is decomposed by Principal Component Analysis (PCA) to find a set of basis eigenvectors such that their linear combination (including mean profile) can result in a frequency-dependent profile template.
        Mean profile ($\tilde{p}_{mean}$) and the basis-eigenvectors ($\hat{e}_i$) are smoothed in the process.
        A template $T(\nu)$ at a particular frequency $\nu$ can be created using the equation:
        \begin{equation}
        T(\nu)=\sum_{i=1}^{n_{eig}}\sum_{j=1}^{n_{B}}c_{ij}B_{tk,j}(\nu)\hat{e}_i+\tilde{p}_{mean}\end{equation}
        where, the first sum runs over the number of basis-eigenvectors $n_{eig}$ and the second over $n_B$, the number of basis splines used in the fit. $\sum_{j=1}^{n_{B}}c_{ij}B_{tk,j}(\nu)$ are the coefficients of the eigenvectors which can evolve with frequency to capture the profile evolution. The default constraint in the software, to limit the number of eigenvectors, is to set a threshold value for S/N of eigenvectors. To determine the threshold value, we first allowed all possible eigenvectors with positive S/N. Then we check each eigenvector's shape and profile evolution of coefficient of eigenvector with frequency. For most eigenvectors, which capture noise, have low S/N with their shapes seem to be noise-like, and, more importantly, their coefficient shows a random behavior with frequency. It leads us to decide S/N cut-offs to determine a reduced set of orthogonal basis-eigenvectors. Note that,  alignment error in the data sets can lead to smearing in the average portrait resulting in more number of eigenvectors. To avoid misalignment we use only clean data sets (with high S/N) to create average portrait and select the initial guess as the highest S/N epoch data. For most of the pulsars we find that only 1-2 components were fitting the above criteria for basis-eigenvectors. Along with the WB ToAs for PTA MSPs, a sample of wide-band timing analysis ``jupyter notebooks'' developed for the uGMRT band-3, band-4, and band-5, which needs a folded data cube in PSRFITS format, are available in the github\footnote{\url{https://github.com/Shyamss6027557/Wide-band-timing-at-Low-frequencies-with-uGMRT}}. Over the uGMRT frequency bands, the pulse profiles show a clear difference in profile shapes from one band to another. So the absolute DM values are expected to be different for the two non-simultaneous bands. Also, it is difficult to completely disentangle DM from frequency evolution of the profile, which can result in different absolute DM values in the NB and WB analysis. To remove these offsets, we've subtracted the weighted mean of all DM values, in a particular band (and particular analysis),
        from all the individual DM values of that band.
 
\section{Results}
\label{sec:Results}

\subsection{Validation of the WB timing pipeline} 
With the aim to validate the WB timing pipeline we have carried out a comparative study of NB and WB timing analysis for the PTA pulsars to obtain DM variations over nine months of timing campaign using band-3 (300-500 MHz) and band-5 (1060-1460 MHz) of the uGMRT. 
 
The pulse profile of PSR J2145$-$0750 has two main components, and the peak amplitude ratio evolves  with frequency (as seen in Fig. \ref{a} and \ref{e}). The pulsar is significantly brighter in band-3 than in band-5, making it one of the best PTA MSPs to follow-up at the low frequencies. We quantify the profile evolution within a frequency band by the number of eigenvectors required to create a frequency-dependent profile template (Section \ref{sec:techniques}).
For PSR J2145$-$0750, 1 and 2 eigenvectors was required to capture most of the profile evolution with frequency in average portrait of band-5 and band-3, respectively.
\begin{figure}[H]
\centering
        \includegraphics[width=\linewidth,keepaspectratio]{J2145-0750_updated.png}
        \vspace{-10mm}
    \caption{Figure shows the DM variation with time for the PTA pulsar J2145$-$0750 in two separate bands of the uGMRT. Top panel: Dark blue and dark red triangles show DM values obtained by NB and WB analysis in band-3. Bottom panel: Top panel points are over-plotted with band-5 DM values. We've used violet and light salmon triangles for NB and WB analysis DM values in band-5. Up triangles are used for coherently de-dispersed data (in band-3) while down triangles for incoherently de-dispersed data (in band-5). Zoomed inset plots for a single epoch DM measurement shows the difference in size of error bars obtained in the NB and WB analysis.}
    \label{DM_curve_J2145-0750}
\end{figure}

For most of the observing epochs in band-3 for this pulsar, the ToA uncertainties from NB and WB analysis are less than 1 $\mu s$. Fig. 2 of \cite{2021ApJS..252....5A} present a median DM and ToA precision, for PSR J2145$-$0750, as 3$\times10^{-4}$\,pc cm$^{-3}\,$ and 200 ns at 820 MHz (with a bandwidth of 186 MHz) with WB analysis using data from the Green Bank Telescope (GBT)\footnote{\url{https://greenbankobservatory.org/science/telescopes/gbt/}}. We estimate the median DM and ToA accuracy of 4$\times10^{-5}$\,pc cm$^{-3}\,$ and 340 ns with WB analysis at band-3 of the uGMRT. The lower frequency pulses arrive at a larger separation which allows more accurate DM measurements at low-frequencies as compared to high radio frequencies. Also, the DM precision is influenced by the steep spectral index of PSR J2145$-$0750 ($\sim$ $-$2 as seen in Table \ref{table_flux}). Unlike \cite{2021ApJS..252....5A}, we have not carried out any noise analysis for modeling red and white noise components of the timing residuals to refine the ToA measurements, which can contribute to the difference of parameter precision achieved in these two works.  Recently, with the uGMRT band-3, \cite{2021arXiv211206908N} reported similar (to this work) DM and ToA precision (after accounting for twice more integration time given per epoch in \cite{2021arXiv211206908N}) over $\sim$ 150 days.  \cite{2021arXiv210105334K} reported inter-band median DM precision of 3$\times10^{-4}$\,pc cm$^{-3}\,$ using simultaneous 400$-$500 MHz (with 16\% of the GMRT array) and 1360$-$1460 MHz (with 40\% of the  GMRT array) observations. We have not re-scaled the parameter uncertainties by the value of minimised chi-square by employing the ``norescale" option in TEMPO2 while fitting DM. In addition to different observational frequency configuration, profile template, this can also influence the DM precision  achieved in our work compared to  \cite{2021arXiv210105334K}.
We find the mean DM values from WB and NB timing analysis as 9.00180 pc cm$^{-3}\,$ and 9.00174 pc cm$^{-3}\,$ respectively, only differing by $\sim$ 1$\sigma$ median DM uncertainties of WB analysis.
Fig. \ref{DM_curve_J2145-0750} presents the temporal variations of DM for PSR J2145$-$0750, obtained with NB and WB analysis, in band-3 and band-5 of the uGMRT. Zoomed inset plots show a visualization of the achieved DM precision. The WB timing analysis, in general, provides a more precise DM compared to the NB timing. We find the same DM variations for NB and WB analysis in both band-3 and band-5, separately. Also, the DM variation obtained in band-3 ($\sim$ 4$\times10^{-5}$\,cm$^{-3}\,$pc) lies within the error bars of band-5 ($\sim$ 9$\times10^{-3}$\,cm$^{-3}\,$pc). Overall, we find a good agreement in the output DM variations from the two analyses, in different uGMRT bands, for PSR J2145$-$0750.  

Unlike PSR J2145$-$0750, PSR J1713$+$0747 has a lower flux density in band-3 compared to band-5 (Table \ref{table_flux}). The pulse profiles in both the bands have a single component with significant pulse broadening due to scattering seen in band-3 (Fig \ref{d}) compared to band-5 (Fig \ref{h}). The pulse profile within band-3 also evolves considerably with frequency due to scattering. We required three eigenvectors for modeling its average portrait in band-3 and one in band-5.
For PSR J1713$+$0747, \cite{2021ApJS..252....5A} reported a median DM and ToA precision of $\sim 1\times 10^{-4}$\,pc cm$^{-3}$ and 90 ns respectively, at 820 MHz  with WB analysis. 
From WB timing analysis of the uGMRT band-3 data, we derive a median DM and ToA precision of 2.4$\times10^{-4}$\,pc cm$^{-3}\,$ and 960 ns, respectively. The NB timing analysis provides lower DM and ToA precision than WB (as seen in Table \ref{Table_2}) with same variation as seen for PSR J2145$-$0750. In addition to pulse broadening, the lower flux density value for PSR J1713$+$0750 at 400 MHz, compared to higher frequencies ({ATNF pulsar catalogue$^{\ref{footnotemark1}}$}), could be limiting the DM precision at 400 MHz.
Fig \ref{b} and \ref{f} show folded pulse profiles for PSR J1640$+$2224 in band-3 and band-5, respectively. The steep spectral nature of this pulsar makes it much brighter in band-3 as compared to band-5. The average portrait in band-3 needs one eigenvector to model its profile evolution with frequency. However, the band-5 average portrait doesn't require any such eigenvector in WB timing analysis.
The WB timing analysis of 820~MHz observations reported in \cite{2021ApJS..252....5A} gives a median DM and ToA uncertainties, for J1640+2224 as $\sim 1\times10^{-4}$\,cm$^{-3}\,$pc and 80 ns, respectively. 
 With the uGMRT band-3 observations, we estimate the median DM and ToA accuracy of 2$\times10^{-5}$\,cm$^{-3}\,$pc and 150 ns respectively, in WB analysis. The steep spectral nature of this pulsar aided with WB timing analysis provides improved DM precision with the uGMRT low-frequency observations.
 PSR J1909$-$3744 has a single component pulse profile both in band-3 (Fig \ref{c}) and band-5 (Fig \ref{g}), and has higher detection significance in band-3 compared to band-5. Three eigenvectors were required to model the profile evolution with frequency in average portrait of band-3. No such eigenvector is required to create its template in band-5.
For PSR J1909$-$3744, the median DM and ToA precision reported by \cite{2021ApJS..252....5A} using 820~MHz observations are $\sim 5\times10^{-5}$\,pc cm$^{-3}\,$ and 30 ns, respectively.
With uGMRT band-3, we estimate the median DM and ToA accuracy of $4\times10^{-5}$\,cm$^{-3}\,$pc and 260 ns using the WB analysis. A timing baseline with longer span, and regular cadence, can provide long-term DM behavior at low frequencies resulting from WB analysis. Moreover, monitoring at band-5 allows comparison with the results with \cite{2021ApJS..252....5A}.  
 \vspace{-3mm}
\subsection{DM variations of GMRT pulsars}
 All four non-PTA pulsars discovered with the GMRT show strong profile evolution with frequency as seen in Fig. \ref{GMRT_pulsar_profiles}. PSR J1646$-$2142 has an interesting profile evolution with frequency. In addition to the frequency evolution of peak amplitude ratio of two components, the component separation in band-3 also decreases at lower frequencies opposite to the radius-to-frequency mapping seen for some pulsars \citep{2012hpa..book.....L}. The peak amplitude ratio of the two components also changes from band-3 to band-4 (Fig. \ref{m}). PSR J2144$-$5237 shows similar evolution of component amplitude with frequency (Fig. \ref{b}). We see a wide-profile single component feature for PSR J1120$-$3618, whereas PSR J1828$+$0625 exhibits a single narrow component pulse profile. We find a single eigenvector to be optimal for profile modeling of average portrait in band-3 for three non-PTA pulsars except J1828+0625, which doesn't require any eigenvector. The band-4 analysis for PSR J2144$-$5237 needs no eigenvector for the average portrait, whereas only one eigenvector is significant in modeling the average portrait of PSR J1646$-$2142. The temporal DM variation of these pulsars, obtained from NB and WB analysis, are shown in Fig. \ref{J1120_GMRT}, \ref{J1646_GMRT}, \ref{J2144_GMRT}, and \ref{J1828_GMRT}.
 The GMRT non-PTA pulsars (J1120$-$3618, J1646$-$2142, J2144$-$5237, and J1828$+$0625) show similar DM variations for NB and WB analysis like PTA MSPs. The median DM uncertainties are around 10$^{-4}$\,pc cm$^{-3}\,$ achieved in band-3 using the WB timing. The mean DM values for these pulsars from WB and NB analysis are:  45.1346 pc cm$^{-3}\,$ and 45.1345 pc cm$^{-3}\,$ for PSR J1120$-$3618; 29.7567 pc cm$^{-3}\,$ and 29.7566 pc cm$^{-3}\,$ for PSR J1646$-$2142; 19.5501 pc cm$^{-3}\,$ for PSR J2144$-$5237; 22.4165 pc cm$^{-3}\,$ and 22.4164 for PSR J1828$+$0625. Thus mean DM offsets between NB and WB timing are insignificant.  
\vspace{-3mm}
\begin{figure}[H]
\centering
        \includegraphics[width=0.9\linewidth,keepaspectratio]{J1120-3618_edited.jpg}
        \vspace{-5mm}
    \caption{DM variation with time for the non-PTA pulsar J1120$-$3618 in band-3 of uGMRT. Dark red triangles show DM values obtained by WB analysis. Dark blue triangles show DM values obtained by NB analysis. Up triangles are used for coherently de-dispersed data while down triangles for incoherently de-dispersed data.}
    \label{J1120_GMRT}
\end{figure}
\vspace{-3mm}
\begin{figure}[H]
\centering
        \includegraphics[width=0.9\linewidth,keepaspectratio]{J1646-2142_edited.png}
                \vspace{-5mm}
    \caption{Figure showing DM variation with time for the non-PTA pulsar J1646$-$2142 in two separate bands of the uGMRT. Dark red and light salmon triangles show DM values obtained by WB analysis in band-3 and band-4, respectively. Dark blue and violet triangles show DM values obtained by NB analysis in band-3 and band-4, respectively. Symbols are same as Fig. 5.}
    \label{J1646_GMRT}
\end{figure}
\vspace{-5mm}
\begin{figure}[H]
\centering
        \includegraphics[width=0.9\linewidth,keepaspectratio]{J2144-5237_edited.png}
        \vspace{-5mm}
    \caption{Figure showing DM variation with time for the non-PTA pulsar J2144$-$5237 in two separate bands of the uGMRT. The colour and symbol schemes are same as in Fig. 6.}
    \label{J2144_GMRT}
\end{figure}
\vspace{-5mm}
\begin{figure}[H]
\centering
        \includegraphics[width=0.9\linewidth,keepaspectratio]{J1828+0625_edited.png}
        \vspace{-5mm}
    \caption{DM variation with time for the non-PTA pulsar J1828$+$0625 in band-3 of the uGMRT. The colour scheme is same as Fig. 5.}
    \label{J1828_GMRT}
\end{figure}
\vspace{-5mm}
As the ToAs and DMs for individual epochs are optimized simultaneously in the WB analysis, they are in general more precisely estimated in WB timing than NB timing for the non-PTA pulsars as well. For PSR J1646$-$2142 and J2144$-$5237, band-4 data gives an order of magnitude less DM precision $\sim$ 10$^{-3}$\,pc cm$^{-3}\,$. Not much DM variation ($<$ $\pm$3$\sigma_{DM}$) is seen for the GMRT pulsar over 2$-$4 years of timing baseline. Only PSR J1120$-$3618 shows $>$ $\pm$4$\sigma_{DM}$ change of DM from 3.25 years of timing follow-up.

\begin{table}[H]
\begin{tabular}{|c |c |c |c |c |c |c |c|}
 \hline
 {{PSR}} & {$\sigma_{ToA}$} & {$\sigma_{ToA}$}& {$timing$} &  {$timing$} & {$\sigma_{DM}$} & {$\sigma_{DM}$} \\ 
 &{{NB}}&{{WB}} & {$residual_{{NB}}$} &  {$residual_{{WB}}$}&{{NB}}&{{WB}} \\ 
  & {$\mu s$} &{$\mu s$}&{$\mu s$}&{$\mu s$}& {$10^{-4}cm^{-3}pc$} & {$10^{-4}cm^{-3}pc$}\\ 
  \hline\hline
  {{J1646-2142}} & {4.29} & {3.13} & {6.4} & {6.1} & {4.5} & {3.6} \\ 
  {{J2144-5237}} & {3.84} & {2.69} & {4.7} & {4.7} & {4.0} & {4.2} \\ 
  {{J1120-3618}} & {7.05}  & {4.71} & {13.2} & {7.0} & {7.3} & {7.7} \\ 
  {{J1828+0625}} & {4.47}  & {4.03} & {10.9} & {7.4} & {5.8} & {6.2} \\
  \hline
  {{J2145-0750}} & {0.67}  & { 0.34} & {-} & {-} & {1.5} & {0.4} \\ 
  {{J1640+2224}} & {0.42}  & {0.15} & {-} & {-} & {0.2} & {0.2} \\ 
  {{J1713+0747}} & {1.95}  & {0.96} & {-} & {-} & {3.2} & {2.4} \\ 
  {{J1909-3744}} & {0.47}  & {0.26} & {-} & {-} & {0.6} & {0.4} \\ 
 \hline
 \end{tabular}
 \caption{The table lists the median ToA, DM precision, and post-fit timing residuals obtained from NB and WB analysis using band-3 data. The available timing span for PTA pulsars was $<1$ year, whereas the non-PTA pulsars were observed for 2$-$4 years.}
 \label{Table_2}
 \end{table}
\vspace{-5mm}
\begin{figure}[H]
    \centering
    \includegraphics[width=0.85\linewidth, keepaspectratio]{Uncertainty_mix.png}
    \caption{Figure showing median ToA (top panel) and DM (bottom panel) uncertainties obtained in the two (NB and WB) analysis for eight pulsars. Error bars represent the range of precision (standard deviation) obtained for the individual pulsar data sets. We used green, blue, and gray colors to represent the values obtained from WB analysis in band-3, 4, and 5, respectively. Similarly, light-green, sky-blue, and light-gray colors are used for NB analysis values in band-3, 4, and 5, respectively. All the PTA pulsars were observed in band-3 and band-5, whereas the non-PTA GMRT pulsars are mostly observed in band-3 (two of the non-PTA pulsars are having good detection in band-4 as well).}
    \label{toa_dm_NB}
\end{figure}

\begin{figure}[H]
    \centering
    \includegraphics[width=0.8\linewidth, keepaspectratio]{Model_parameters_new_1.jpg}
    \caption{Figure showing a comparison of the timing models of 4 non-PTA GMRT MSPs, in band-3, from NB and WB timing analysis. The x-axis shows the differences of the fitted parameter (astrometric, spin, and binary) values normalized by the uncertainties from the NB timing model, i.e., $(X_{NB}-X_{WB})/\sigma_X^{NB}$ where X is the pulsar's model parameter. The error bars have a length equal to the ratio of parameter uncertainties from WB and NB models, i.e., $\sigma_X^{WB}/\sigma_X^{NB}$. Model parameters for all MSPs lie within $\pm 3 \sigma_{NB}$ for the two analyses. This is similar to Fig. 3 of \cite{2021ApJS..252....5A}.}
    \label{model_parameters}
\end{figure}
Table \ref{Table_2} presents a comparison of results from NB and WB timing analysis which lists ToA, DM precision as well as post-fit timing residuals using band-3 observations. Fig. \ref{toa_dm_NB} shows median ToA and DM accuracy along with their standard deviations (plotted as error bars) in NB and WB analysis, respectively, for all pulsars in band-3/4/5.
The WB technique gives better DM and ToA precision than NB for all eight pulsars. We applied WB technique to 2-4 years of GMRT band-3 data to achieve post-fit residuals of $<$ 10 $\mu$s for all non-PTA pulsars. For the PTA pulsars, we have not fitted for long-term timing model (except DM for individual epochs), due to the availability of shorter timing span ($<$ year) with sparse sampling. The ephimeries, for the PTA pulsars, are obtained from the NANOGrav archive\footnote{\url{https://data.nanograv.org/}}. The results from WB timing presented in Fig. \ref{toa_dm_NB} clearly show that band-3 of the uGMRT is the optimal frequency for precision DM measurements. However, the ToA uncertainties are not significantly different in all three frequency bands for most of the observed pulsars. 

In Fig. \ref{model_parameters} we compare the long-term timing models of the four non-PTA pulsars from NB and WB timing analysis. We plot the differences of the fitted parameter values normalized by the uncertainties from the narrow-band timing model ($\sigma_{NB}$). The error bars have a length of $\sigma_{WB}$/$\sigma_{NB}$. We find that the model difference of these two timing methods is within the $\pm$ 3-$\sigma$ of $\sigma_{NB}$ confirming generation of similar long-term timing models from both the techniques (i.e. NB or WB). Moreover, in general, WB timing is providing better precision for the fitted model parameters.


\section{Conclusions}
\label{sec:conclusions}
We provide a comparative study of a wide-band timing analysis with the conventional narrow-band timing analysis at low frequency using the uGMRT for a set of non-PTA GMRT MSPs and for some well studied PTA MSPs. The WB timing provides  $\mu$s timing precision as well as 10$^{-4}$ pc cm$^{-3}$ DM measurements for non-PTA pulsars in 300$-$500 MHz. For PTA pulsars, we achieve sub-$\mu$s timing and 10$^{-5}$ pc cm$^{-3}$ DM precision. The non-PTA pulsars have flux densities $\sim$ mJy, similar to the MSPs expected to be found in ongoing surveys. Thus this work illustrates the typical DM and timing precision that can be achieved for newly discovered pulsars from low-frequency follow-up studies.


The fitted model parameters from NB and WB analyses for the non-PTA pulsars in 2$-$4 years of long-term timing agree well within $\pm$ 3 $\sigma$ uncertainties confirming the applicability of WB analysis for long-term timing. We are getting up to a factor of $\sim$ 2-times better, post-fit timing residuals in WB timing compared to NB analysis. Similarly, WB timing generally gives an improved DM precision for all observed pulsars. The mean DM values obtained from the two methods (NB and WB), for the non-PTA pulsars in band-3, and 4 are within the DM errors.  No significant DM variations observed for the GMRT pulsars, except PSR J1120$-$3618. The DM variation for the PTA MSP presented in this paper J2145$-$0750, also agrees well with \cite{2021arXiv211206908N} having first 150 days of overlapping time span, but with different observing configuration. The mean DM values from NB and WB timing for PSR J2145$-$0750 are within median DM uncertainties. We also report the WB timing results for PTA MSPs in band-5 (1060$-$1460 MHz) and for non-PTA MSPs in band-4 (550$-$750 MHz) as well. Both these bands give DM precision in the order of 10$^{-3}$ pc cm$^{-3}$. Thus band-3 of the GMRT is providing much higher DM precision for all the eight MSPs compared to the other observing bands.      

For the commonly observed PTA pulsars, the DM precision obtained from the 300$-$500 MHz bands are at par with the values values reported in   \cite{2021arXiv211206908N} and in some cases better than \cite{2021ApJS..252....5A}. In the case of the non-PTA GMRT pulsars even with an order of magnitude lower flux densities than the PTA MSPs, the achieved DM precision (Table \ref{Table_2}), is comparable with the higher frequency measurements for PTA MSPs \citep{2021ApJS..252....5A} making them as potential candidates to include in the IPTA experiment in a search for a GW background aided by the more sensitive upcoming telescopes providing better ToA precision. Since at the intermediate signal regime of stochastic background of GWs, the detection significance strongly depends on the number of pulsars included in the array and the timing precision \citep{2013CQGra..30v4015S}, such low frequency follow-up aided with wide-band timing can play an important role.
Following the work by \cite{2021arXiv211206908N} for PTA MSPs and the current work for newly discovered MSPs, the prospect of using a low frequency observing facility at a sensitive telescope like the uGMRT for high precision timing studies to aid the global PTA efforts is clearly evident. 

Moreover, the timing with the full GMRT array (70\%) in band-3 (with a gain of 7 K/Jy), presented here, is complementary to the currently existing wide-band timing facilities like MeerKAT and CHIME provoding lowest frequency coverage of 580$-$1670 MHz \citep{2020PASA...37...28B} and 400$-$800 MHz \citep{2021ApJS..255....5C}, respectively.

\section{Acknowledgments}
We acknowledge the support of the Department of Atomic Energy, Government of India,
under project no. 12-R$\&$D-TFR-5.02-0700. The GMRT is run by the institute National Centre for
Radio Astrophysics of the Tata Institute of Fundamental Research, India. We acknowledge the support of GMRT telescope operators for observations. 

